diff --git a/.travis.yml b/.travis.yml
index e137af3..57d01a4 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -1,7 +1,9 @@
 language: go
 go:
- - 1.7
- - 1.8
+ - "1.7"
+ - "1.8"
+ - "1.9"
+ - "1.10"
 env:
  # Temporary workaround for go 1.6
  - GODEBUG=cgocheck=0
diff --git a/base/classifier.go b/base/classifier.go
index 5f61cf2..11adff9 100644
--- a/base/classifier.go
+++ b/base/classifier.go
@@ -1,7 +1,7 @@
 package base
 
 import (
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 // Classifier implementations predict categorical class labels.
@@ -37,7 +37,7 @@ type BaseClassifier struct {
 }
 
 type BaseRegressor struct {
-	Data   mat64.Dense
+	Data   mat.Dense
 	Name   string
 	Labels []float64
 }
diff --git a/base/conversion.go b/base/conversion.go
index 03fcf62..78d344e 100644
--- a/base/conversion.go
+++ b/base/conversion.go
@@ -3,7 +3,7 @@ package base
 import (
 	"fmt"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 func checkAllAttributesAreFloat(attrs []Attribute) error {
@@ -18,8 +18,8 @@ func checkAllAttributesAreFloat(attrs []Attribute) error {
 
 // ConvertRowToMat64 takes a list of Attributes, a FixedDataGrid
 // and a row number, and returns the float values of that row
-// in a mat64.Dense format.
-func ConvertRowToMat64(attrs []Attribute, f FixedDataGrid, r int) (*mat64.Dense, error) {
+// in a mat.Dense format.
+func ConvertRowToMat64(attrs []Attribute, f FixedDataGrid, r int) (*mat.Dense, error) {
 
 	err := checkAllAttributesAreFloat(attrs)
 	if err != nil {
@@ -27,7 +27,7 @@ func ConvertRowToMat64(attrs []Attribute, f FixedDataGrid, r int) (*mat64.Dense,
 	}
 
 	// Allocate the return value
-	ret := mat64.NewDense(1, len(attrs), nil)
+	ret := mat.NewDense(1, len(attrs), nil)
 
 	// Resolve all the attributes
 	attrSpecs := ResolveAttributes(f, attrs)
@@ -42,8 +42,8 @@ func ConvertRowToMat64(attrs []Attribute, f FixedDataGrid, r int) (*mat64.Dense,
 }
 
 // ConvertAllRowsToMat64 takes a list of Attributes and returns a vector
-// of all rows in a mat64.Dense format.
-func ConvertAllRowsToMat64(attrs []Attribute, f FixedDataGrid) ([]*mat64.Dense, error) {
+// of all rows in a mat.Dense format.
+func ConvertAllRowsToMat64(attrs []Attribute, f FixedDataGrid) ([]*mat.Dense, error) {
 
 	// Check for floats
 	err := checkAllAttributesAreFloat(attrs)
@@ -53,14 +53,14 @@ func ConvertAllRowsToMat64(attrs []Attribute, f FixedDataGrid) ([]*mat64.Dense,
 
 	// Return value
 	_, rows := f.Size()
-	ret := make([]*mat64.Dense, rows)
+	ret := make([]*mat.Dense, rows)
 
 	// Resolve all attributes
 	attrSpecs := ResolveAttributes(f, attrs)
 
 	// Set the values in each return value
 	for i := 0; i < rows; i++ {
-		cur := mat64.NewDense(1, len(attrs), nil)
+		cur := mat.NewDense(1, len(attrs), nil)
 		for j, a := range attrSpecs {
 			cur.Set(0, j, UnpackBytesToFloat(f.Get(a, i)))
 		}
diff --git a/base/domain.go b/base/domain.go
index 224aece..357efc5 100644
--- a/base/domain.go
+++ b/base/domain.go
@@ -7,7 +7,7 @@ import (
 	"encoding/gob"
 	"io/ioutil"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 // An Estimator is object that can ingest some data and train on it.
@@ -27,7 +27,7 @@ type Model interface {
 }
 
 type BaseEstimator struct {
-	Data *mat64.Dense
+	Data *mat.Dense
 }
 
 // SaveEstimatorToGob serialises an estimator to a provided filepath, in gob format.
diff --git a/base/mat.go b/base/mat.go
index 02beb9e..f644cf3 100644
--- a/base/mat.go
+++ b/base/mat.go
@@ -3,18 +3,18 @@ package base
 import (
 	"bytes"
 	"fmt"
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type Mat64Instances struct {
 	attributes []Attribute
 	classAttrs map[int]bool
-	Data       *mat64.Dense
+	Data       *mat.Dense
 	rows       int
 }
 
 // InstancesFromMat64 returns a new Mat64Instances from a literal provided.
-func InstancesFromMat64(rows, cols int, data *mat64.Dense) *Mat64Instances {
+func InstancesFromMat64(rows, cols int, data *mat.Dense) *Mat64Instances {
 
 	var ret Mat64Instances
 	for i := 0; i < cols; i++ {
diff --git a/base/mat_test.go b/base/mat_test.go
index f4b8de2..b2ed083 100644
--- a/base/mat_test.go
+++ b/base/mat_test.go
@@ -1,15 +1,15 @@
 package base
 
 import (
-	"github.com/gonum/matrix/mat64"
 	. "github.com/smartystreets/goconvey/convey"
+	"gonum.org/v1/gonum/mat"
 	"testing"
 )
 
 func TestInlineMat64Creation(t *testing.T) {
 
 	Convey("Given a literal array...", t, func() {
-		mat := mat64.NewDense(4, 3, []float64{
+		mat := mat.NewDense(4, 3, []float64{
 			1, 0, 1,
 			0, 1, 1,
 			0, 0, 0,
diff --git a/base/serialize.go b/base/serialize.go
index f1e8046..7f4ba7b 100644
--- a/base/serialize.go
+++ b/base/serialize.go
@@ -72,7 +72,7 @@ func tarPrefix(prefix string, suffix string) string {
 	if prefix == "" {
 		return suffix
 	}
-	return fmt.Sprintf("%s/%s")
+	return fmt.Sprintf("%s/%s", prefix, suffix)
 }
 
 // ClassifierMetadataV1 is what gets written into METADATA
diff --git a/base/serialize_instances.go b/base/serialize_instances.go
index 278f01e..7239776 100644
--- a/base/serialize_instances.go
+++ b/base/serialize_instances.go
@@ -164,16 +164,19 @@ func DeserializeInstancesFromTarReader(tr *FunctionalTarReader, prefix string) (
 
 	// Finally, read the values out of the data section
 	for i := 0; i < rowCount; i++ {
-		for _, s := range specs {
+		for j, s := range specs {
 			r := ret.Get(s, i)
 			n, err := reader.Read(r)
 			if n != len(r) {
-				return nil, fmt.Errorf("Expected %d bytes (read %d) on row %d", len(r), n, i)
-			}
-			if err != nil {
-				return nil, fmt.Errorf("Read error: %s", err)
+				return nil, WrapError(fmt.Errorf("Expected %d bytes (read %d) on row %d", len(r), n, i))
 			}
 			ret.Set(s, i, r)
+			if err != nil {
+				if i == rowCount-1 && j == len(specs)-1 && err == io.EOF {
+					break
+				}
+				return nil, WrapError(fmt.Errorf("Read error in data section (at row %d from %d, attr %d from %d): %s", i, rowCount, j, len(specs), err))
+			}
 		}
 	}
 
@@ -300,6 +303,9 @@ func SerializeInstancesToTarWriter(inst FixedDataGrid, tw *tar.Writer, prefix st
 	}
 
 	allSpecs := ResolveAttributes(inst, allAttrs)
+	if len(allSpecs) != len(allAttrs) {
+		return WrapError(fmt.Errorf("Error resolving all Attributes: resolved %d, expected %d", len(allSpecs), len(allAttrs)))
+	}
 
 	// First, estimate the amount of data we'll need...
 	dataLength := int64(0)
diff --git a/clustering/dbscan.go b/clustering/dbscan.go
index a40ee79..67a752e 100644
--- a/clustering/dbscan.go
+++ b/clustering/dbscan.go
@@ -1,7 +1,7 @@
 package clustering
 
 import (
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	"github.com/sjwhitworth/golearn/base"
 	"github.com/sjwhitworth/golearn/metrics/pairwise"
 	"math/big"
@@ -22,7 +22,7 @@ type DBSCANParameters struct {
 	MinCount int
 }
 
-func regionQuery(p int, ret *big.Int, dist *mat64.Dense, eps float64) *big.Int {
+func regionQuery(p int, ret *big.Int, dist *mat.Dense, eps float64) *big.Int {
 	rows, _ := dist.Dims()
 	// Return any points within the Eps neighbourhood
 	for i := 0; i < rows; i++ {
@@ -33,7 +33,7 @@ func regionQuery(p int, ret *big.Int, dist *mat64.Dense, eps float64) *big.Int {
 	return ret
 }
 
-func computePairwiseDistances(inst base.FixedDataGrid, attrs []base.Attribute, metric pairwise.PairwiseDistanceFunc) (*mat64.Dense, error) {
+func computePairwiseDistances(inst base.FixedDataGrid, attrs []base.Attribute, metric pairwise.PairwiseDistanceFunc) (*mat.Dense, error) {
 	// Compute pair-wise distances
 	// First convert everything to floats
 	mats, err := base.ConvertAllRowsToMat64(attrs, inst)
@@ -43,7 +43,7 @@ func computePairwiseDistances(inst base.FixedDataGrid, attrs []base.Attribute, m
 
 	// Next, do an n^2 computation of all pairwise distances
 	_, rows := inst.Size()
-	dist := mat64.NewDense(rows, rows, nil)
+	dist := mat.NewDense(rows, rows, nil)
 	for i := 0; i < rows; i++ {
 		for j := i + 1; j < rows; j++ {
 			d := metric.Distance(mats[i], mats[j])
diff --git a/clustering/dbscan_test.go b/clustering/dbscan_test.go
index 182fca3..b2c7eb1 100644
--- a/clustering/dbscan_test.go
+++ b/clustering/dbscan_test.go
@@ -2,7 +2,7 @@ package clustering
 
 import (
 	"bufio"
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	"github.com/sjwhitworth/golearn/base"
 	"github.com/sjwhitworth/golearn/metrics/pairwise"
 	. "github.com/smartystreets/goconvey/convey"
@@ -76,8 +76,8 @@ func TestDBSCANDistanceMetric(t *testing.T) {
 
 	Convey("Check the distance function is sane...", t, func() {
 
-		d1 := mat64.NewDense(1, 2, nil)
-		d2 := mat64.NewDense(1, 2, nil)
+		d1 := mat.NewDense(1, 2, nil)
+		d2 := mat.NewDense(1, 2, nil)
 
 		d1.Set(0, 0, 0.494260967249)
 		d1.Set(0, 1, 1.45106696541)
diff --git a/kdtree/kdtree.go b/kdtree/kdtree.go
index e828d69..d39ed08 100644
--- a/kdtree/kdtree.go
+++ b/kdtree/kdtree.go
@@ -2,7 +2,7 @@ package kdtree
 
 import (
 	"errors"
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	"github.com/sjwhitworth/golearn/metrics/pairwise"
 	"sort"
 )
@@ -140,8 +140,8 @@ func (t *Tree) Search(k int, disType pairwise.PairwiseDistanceFunc, target []flo
 
 func (t *Tree) searchHandle(k int, disType pairwise.PairwiseDistanceFunc, target []float64, h *heap, n *node) {
 	if n.feature == -1 {
-		vectorX := mat64.NewDense(len(target), 1, target)
-		vectorY := mat64.NewDense(len(target), 1, n.value)
+		vectorX := mat.NewDense(len(target), 1, target)
+		vectorY := mat.NewDense(len(target), 1, n.value)
 		length := disType.Distance(vectorX, vectorY)
 		h.insert(n.value, length, n.srcRowNo)
 		return
@@ -157,8 +157,8 @@ func (t *Tree) searchHandle(k int, disType pairwise.PairwiseDistanceFunc, target
 		t.searchHandle(k, disType, target, h, n.right)
 	}
 
-	vectorX := mat64.NewDense(len(target), 1, target)
-	vectorY := mat64.NewDense(len(target), 1, n.value)
+	vectorX := mat.NewDense(len(target), 1, target)
+	vectorY := mat.NewDense(len(target), 1, n.value)
 	length := disType.Distance(vectorX, vectorY)
 
 	if k > h.size() {
@@ -177,8 +177,8 @@ func (t *Tree) searchHandle(k int, disType pairwise.PairwiseDistanceFunc, target
 			t.searchAllNodes(k, disType, target, h, n.left)
 		}
 	} else {
-		vectorX = mat64.NewDense(1, 1, []float64{target[n.feature]})
-		vectorY = mat64.NewDense(1, 1, []float64{n.value[n.feature]})
+		vectorX = mat.NewDense(1, 1, []float64{target[n.feature]})
+		vectorY = mat.NewDense(1, 1, []float64{n.value[n.feature]})
 		length = disType.Distance(vectorX, vectorY)
 
 		if h.maximum().length > length {
@@ -192,8 +192,8 @@ func (t *Tree) searchHandle(k int, disType pairwise.PairwiseDistanceFunc, target
 }
 
 func (t *Tree) searchAllNodes(k int, disType pairwise.PairwiseDistanceFunc, target []float64, h *heap, n *node) {
-	vectorX := mat64.NewDense(len(target), 1, target)
-	vectorY := mat64.NewDense(len(target), 1, n.value)
+	vectorX := mat.NewDense(len(target), 1, target)
+	vectorY := mat.NewDense(len(target), 1, n.value)
 	length := disType.Distance(vectorX, vectorY)
 
 	if k > h.size() {
diff --git a/knn/knn.go b/knn/knn.go
index 5554863..c8fa9cc 100644
--- a/knn/knn.go
+++ b/knn/knn.go
@@ -8,7 +8,7 @@ import (
 	"fmt"
 
 	"github.com/gonum/matrix"
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	"github.com/sjwhitworth/golearn/base"
 	"github.com/sjwhitworth/golearn/kdtree"
 	"github.com/sjwhitworth/golearn/metrics/pairwise"
@@ -434,11 +434,11 @@ func (KNN *KNNRegressor) Fit(values []float64, numbers []float64, rows int, cols
 		panic(matrix.ErrShape)
 	}
 
-	KNN.Data = mat64.NewDense(rows, cols, numbers)
+	KNN.Data = mat.NewDense(rows, cols, numbers)
 	KNN.Values = values
 }
 
-func (KNN *KNNRegressor) Predict(vector *mat64.Dense, K int) float64 {
+func (KNN *KNNRegressor) Predict(vector *mat.Dense, K int) float64 {
 	// Get the number of rows
 	rows, _ := KNN.Data.Dims()
 	rownumbers := make(map[int]float64)
diff --git a/linear_models/linear_regression.go b/linear_models/linear_regression.go
index 608a87f..b726539 100644
--- a/linear_models/linear_regression.go
+++ b/linear_models/linear_regression.go
@@ -7,7 +7,7 @@ import (
 
 	"fmt"
 	_ "github.com/gonum/blas"
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 var (
@@ -59,8 +59,8 @@ func (lr *LinearRegression) Fit(inst base.FixedDataGrid) error {
 
 	// Split into two matrices, observed results (dependent variable y)
 	// and the explanatory variables (X) - see http://en.wikipedia.org/wiki/Linear_regression
-	observed := mat64.NewDense(rows, 1, nil)
-	explVariables := mat64.NewDense(rows, cols, nil)
+	observed := mat.NewDense(rows, 1, nil)
+	explVariables := mat.NewDense(rows, cols, nil)
 
 	// Build the observed matrix
 	inst.MapOverRows(classAttrSpecs, func(row [][]byte, i int) (bool, error) {
@@ -80,13 +80,13 @@ func (lr *LinearRegression) Fit(inst base.FixedDataGrid) error {
 	})
 
 	n := cols
-	qr := new(mat64.QR)
+	qr := new(mat.QR)
 	qr.Factorize(explVariables)
-	var q, reg mat64.Dense
-	q.QFromQR(qr)
-	reg.RFromQR(qr)
+	var q, reg mat.Dense
+	qr.QTo(&q)
+	qr.RTo(&reg)
 
-	var transposed, qty mat64.Dense
+	var transposed, qty mat.Dense
 	transposed.Clone(q.T())
 	qty.Mul(&transposed, observed)
 
diff --git a/meta/one_v_all.go b/meta/one_v_all.go
index ded4239..e6a7e6f 100644
--- a/meta/one_v_all.go
+++ b/meta/one_v_all.go
@@ -176,7 +176,7 @@ func (m *OneVsAllModel) LoadWithPrefix(reader *base.ClassifierDeserializer, pref
 		attrMap := make(map[base.Attribute]base.Attribute)
 
 		for j := 0; j < int(numAttrsInMapU64); j++ {
-			mapTupleKey := reader.Prefix(mapPrefix, fmt.Sprintf("%d"))
+			mapTupleKey := reader.Prefix(mapPrefix, fmt.Sprintf("%d", j))
 			mapKeyKeyKey := reader.Prefix(mapTupleKey, "KEY")
 			mapKeyValKey := reader.Prefix(mapTupleKey, "VAL")
 
@@ -289,7 +289,7 @@ func (m *OneVsAllModel) SaveWithPrefix(writer *base.ClassifierSerializer, prefix
 		}
 		j := 0
 		for key := range f.attrs {
-			mapTupleKey := writer.Prefix(mapPrefix, fmt.Sprintf("%d"))
+			mapTupleKey := writer.Prefix(mapPrefix, fmt.Sprintf("%d", j))
 			mapKeyKeyKey := writer.Prefix(mapTupleKey, "KEY")
 			mapKeyValKey := writer.Prefix(mapTupleKey, "VAL")
 
diff --git a/metrics/pairwise/chebyshev.go b/metrics/pairwise/chebyshev.go
index c169e83..23a6fb1 100644
--- a/metrics/pairwise/chebyshev.go
+++ b/metrics/pairwise/chebyshev.go
@@ -4,7 +4,7 @@ import (
 	"math"
 
 	"github.com/gonum/matrix"
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type Chebyshev struct{}
@@ -13,7 +13,7 @@ func NewChebyshev() *Chebyshev {
 	return &Chebyshev{}
 }
 
-func (c *Chebyshev) Distance(vectorX *mat64.Dense, vectorY *mat64.Dense) float64 {
+func (c *Chebyshev) Distance(vectorX *mat.Dense, vectorY *mat.Dense) float64 {
 	r1, c1 := vectorX.Dims()
 	r2, c2 := vectorY.Dims()
 	if r1 != r2 || c1 != c2 {
diff --git a/metrics/pairwise/chebyshev_test.go b/metrics/pairwise/chebyshev_test.go
index 61b6a84..2fe39fc 100644
--- a/metrics/pairwise/chebyshev_test.go
+++ b/metrics/pairwise/chebyshev_test.go
@@ -3,17 +3,17 @@ package pairwise
 import (
 	"testing"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	. "github.com/smartystreets/goconvey/convey"
 )
 
 func TestChebyshev(t *testing.T) {
-	var vectorX, vectorY *mat64.Dense
+	var vectorX, vectorY *mat.Dense
 	chebyshev := NewChebyshev()
 
 	Convey("Given two vectors", t, func() {
-		vectorX = mat64.NewDense(4, 1, []float64{1, 2, 3, 4})
-		vectorY = mat64.NewDense(4, 1, []float64{-5, -6, 7, 8})
+		vectorX = mat.NewDense(4, 1, []float64{1, 2, 3, 4})
+		vectorY = mat.NewDense(4, 1, []float64{-5, -6, 7, 8})
 
 		Convey("When calculating distance with two vectors", func() {
 			result := chebyshev.Distance(vectorX, vectorY)
diff --git a/metrics/pairwise/cosine.go b/metrics/pairwise/cosine.go
index 126d749..c587150 100644
--- a/metrics/pairwise/cosine.go
+++ b/metrics/pairwise/cosine.go
@@ -3,7 +3,7 @@ package pairwise
 import (
 	"math"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type Cosine struct{}
@@ -13,17 +13,17 @@ func NewCosine() *Cosine {
 }
 
 // Dot computes dot value of vectorX and vectorY.
-func (c *Cosine) Dot(vectorX *mat64.Dense, vectorY *mat64.Dense) float64 {
-	subVector := mat64.NewDense(0, 0, nil)
+func (c *Cosine) Dot(vectorX *mat.Dense, vectorY *mat.Dense) float64 {
+	subVector := mat.NewDense(0, 0, nil)
 	subVector.MulElem(vectorX, vectorY)
-	result := mat64.Sum(subVector)
+	result := mat.Sum(subVector)
 
 	return result
 }
 
 // Distance computes Cosine distance.
 // It will return distance which represented as 1-cos() (ranged from 0 to 2).
-func (c *Cosine) Distance(vectorX *mat64.Dense, vectorY *mat64.Dense) float64 {
+func (c *Cosine) Distance(vectorX *mat.Dense, vectorY *mat.Dense) float64 {
 	dotXY := c.Dot(vectorX, vectorY)
 	lengthX := math.Sqrt(c.Dot(vectorX, vectorX))
 	lengthY := math.Sqrt(c.Dot(vectorY, vectorY))
diff --git a/metrics/pairwise/cosine_test.go b/metrics/pairwise/cosine_test.go
index 2a3078e..ed74b31 100644
--- a/metrics/pairwise/cosine_test.go
+++ b/metrics/pairwise/cosine_test.go
@@ -3,17 +3,17 @@ package pairwise
 import (
 	"testing"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	. "github.com/smartystreets/goconvey/convey"
 )
 
 func TestCosine(t *testing.T) {
-	var vectorX, vectorY *mat64.Dense
+	var vectorX, vectorY *mat.Dense
 	cosine := NewCosine()
 
 	Convey("Given two vectors", t, func() {
-		vectorX = mat64.NewDense(3, 1, []float64{1, 2, 3})
-		vectorY = mat64.NewDense(3, 1, []float64{2, 4, 6})
+		vectorX = mat.NewDense(3, 1, []float64{1, 2, 3})
+		vectorY = mat.NewDense(3, 1, []float64{2, 4, 6})
 
 		Convey("When doing inner Dot", func() {
 			result := cosine.Dot(vectorX, vectorY)
diff --git a/metrics/pairwise/cranberra.go b/metrics/pairwise/cranberra.go
index 06a99b4..eb36f6c 100644
--- a/metrics/pairwise/cranberra.go
+++ b/metrics/pairwise/cranberra.go
@@ -4,7 +4,7 @@ import (
 	"math"
 
 	"github.com/gonum/matrix"
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type Cranberra struct{}
@@ -20,7 +20,7 @@ func cranberraDistanceStep(num float64, denom float64) float64 {
 	return num / denom
 }
 
-func (c *Cranberra) Distance(vectorX *mat64.Dense, vectorY *mat64.Dense) float64 {
+func (c *Cranberra) Distance(vectorX *mat.Dense, vectorY *mat.Dense) float64 {
 	r1, c1 := vectorX.Dims()
 	r2, c2 := vectorY.Dims()
 	if r1 != r2 || c1 != c2 {
diff --git a/metrics/pairwise/cranberra_test.go b/metrics/pairwise/cranberra_test.go
index fe9a000..147fc25 100644
--- a/metrics/pairwise/cranberra_test.go
+++ b/metrics/pairwise/cranberra_test.go
@@ -3,16 +3,16 @@ package pairwise
 import (
 	"testing"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	. "github.com/smartystreets/goconvey/convey"
 )
 
 func TestCranberrra(t *testing.T) {
-	var vectorX, vectorY *mat64.Dense
+	var vectorX, vectorY *mat.Dense
 	cranberra := NewCranberra()
 
 	Convey("Given two vectors that are same", t, func() {
-		vec := mat64.NewDense(7, 1, []float64{0, 1, -2, 3.4, 5, -6.7, 89})
+		vec := mat.NewDense(7, 1, []float64{0, 1, -2, 3.4, 5, -6.7, 89})
 		distance := cranberra.Distance(vec, vec)
 
 		Convey("The result should be 0", func() {
@@ -21,8 +21,8 @@ func TestCranberrra(t *testing.T) {
 	})
 
 	Convey("Given two vectors", t, func() {
-		vectorX = mat64.NewDense(5, 1, []float64{1, 2, 3, 4, 9})
-		vectorY = mat64.NewDense(5, 1, []float64{-5, -6, 7, 4, 3})
+		vectorX = mat.NewDense(5, 1, []float64{1, 2, 3, 4, 9})
+		vectorY = mat.NewDense(5, 1, []float64{-5, -6, 7, 4, 3})
 
 		Convey("When calculating distance with two vectors", func() {
 			result := cranberra.Distance(vectorX, vectorY)
diff --git a/metrics/pairwise/euclidean.go b/metrics/pairwise/euclidean.go
index 999e2a2..a95cdba 100644
--- a/metrics/pairwise/euclidean.go
+++ b/metrics/pairwise/euclidean.go
@@ -3,7 +3,7 @@ package pairwise
 import (
 	"math"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type Euclidean struct{}
@@ -13,17 +13,17 @@ func NewEuclidean() *Euclidean {
 }
 
 // InnerProduct computes a Eucledian inner product.
-func (e *Euclidean) InnerProduct(vectorX *mat64.Dense, vectorY *mat64.Dense) float64 {
-	subVector := mat64.NewDense(0, 0, nil)
+func (e *Euclidean) InnerProduct(vectorX *mat.Dense, vectorY *mat.Dense) float64 {
+	subVector := mat.NewDense(0, 0, nil)
 	subVector.MulElem(vectorX, vectorY)
-	result := mat64.Sum(subVector)
+	result := mat.Sum(subVector)
 
 	return result
 }
 
 // Distance computes Euclidean distance (also known as L2 distance).
-func (e *Euclidean) Distance(vectorX *mat64.Dense, vectorY *mat64.Dense) float64 {
-	subVector := mat64.NewDense(0, 0, nil)
+func (e *Euclidean) Distance(vectorX *mat.Dense, vectorY *mat.Dense) float64 {
+	subVector := mat.NewDense(0, 0, nil)
 	subVector.Sub(vectorX, vectorY)
 
 	result := e.InnerProduct(subVector, subVector)
diff --git a/metrics/pairwise/euclidean_test.go b/metrics/pairwise/euclidean_test.go
index d022ae6..0996d5a 100644
--- a/metrics/pairwise/euclidean_test.go
+++ b/metrics/pairwise/euclidean_test.go
@@ -3,17 +3,17 @@ package pairwise
 import (
 	"testing"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	. "github.com/smartystreets/goconvey/convey"
 )
 
 func TestEuclidean(t *testing.T) {
-	var vectorX, vectorY *mat64.Dense
+	var vectorX, vectorY *mat.Dense
 	euclidean := NewEuclidean()
 
 	Convey("Given two vectors", t, func() {
-		vectorX = mat64.NewDense(3, 1, []float64{1, 2, 3})
-		vectorY = mat64.NewDense(3, 1, []float64{2, 4, 5})
+		vectorX = mat.NewDense(3, 1, []float64{1, 2, 3})
+		vectorY = mat.NewDense(3, 1, []float64{2, 4, 5})
 
 		Convey("When doing inner product", func() {
 			result := euclidean.InnerProduct(vectorX, vectorY)
diff --git a/metrics/pairwise/manhattan.go b/metrics/pairwise/manhattan.go
index c9d3b75..29c9ba9 100644
--- a/metrics/pairwise/manhattan.go
+++ b/metrics/pairwise/manhattan.go
@@ -4,7 +4,7 @@ import (
 	"math"
 
 	"github.com/gonum/matrix"
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type Manhattan struct{}
@@ -15,7 +15,7 @@ func NewManhattan() *Manhattan {
 
 // Distance computes the Manhattan distance, also known as L1 distance.
 // == the sum of the absolute values of elements.
-func (m *Manhattan) Distance(vectorX *mat64.Dense, vectorY *mat64.Dense) float64 {
+func (m *Manhattan) Distance(vectorX *mat.Dense, vectorY *mat.Dense) float64 {
 	r1, c1 := vectorX.Dims()
 	r2, c2 := vectorY.Dims()
 	if r1 != r2 || c1 != c2 {
diff --git a/metrics/pairwise/manhattan_test.go b/metrics/pairwise/manhattan_test.go
index bd7feba..ea2e6f7 100644
--- a/metrics/pairwise/manhattan_test.go
+++ b/metrics/pairwise/manhattan_test.go
@@ -3,16 +3,16 @@ package pairwise
 import (
 	"testing"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	. "github.com/smartystreets/goconvey/convey"
 )
 
 func TestManhattan(t *testing.T) {
-	var vectorX, vectorY *mat64.Dense
+	var vectorX, vectorY *mat.Dense
 	manhattan := NewManhattan()
 
 	Convey("Given two vectors that are same", t, func() {
-		vec := mat64.NewDense(7, 1, []float64{0, 1, -2, 3.4, 5, -6.7, 89})
+		vec := mat.NewDense(7, 1, []float64{0, 1, -2, 3.4, 5, -6.7, 89})
 		distance := manhattan.Distance(vec, vec)
 
 		Convey("The result should be 0", func() {
@@ -21,8 +21,8 @@ func TestManhattan(t *testing.T) {
 	})
 
 	Convey("Given two vectors", t, func() {
-		vectorX = mat64.NewDense(3, 1, []float64{2, 2, 3})
-		vectorY = mat64.NewDense(3, 1, []float64{1, 4, 5})
+		vectorX = mat.NewDense(3, 1, []float64{2, 2, 3})
+		vectorY = mat.NewDense(3, 1, []float64{1, 4, 5})
 
 		Convey("When calculating distance with column vectors", func() {
 			result := manhattan.Distance(vectorX, vectorY)
diff --git a/metrics/pairwise/pairwise.go b/metrics/pairwise/pairwise.go
index 3a88e67..4dcd453 100644
--- a/metrics/pairwise/pairwise.go
+++ b/metrics/pairwise/pairwise.go
@@ -2,9 +2,9 @@
 package pairwise
 
 import (
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type PairwiseDistanceFunc interface {
-	Distance(vectorX *mat64.Dense, vectorY *mat64.Dense) float64
+	Distance(vectorX *mat.Dense, vectorY *mat.Dense) float64
 }
diff --git a/metrics/pairwise/poly_kernel.go b/metrics/pairwise/poly_kernel.go
index 709742f..ba6d453 100644
--- a/metrics/pairwise/poly_kernel.go
+++ b/metrics/pairwise/poly_kernel.go
@@ -3,7 +3,7 @@ package pairwise
 import (
 	"math"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type PolyKernel struct {
@@ -17,18 +17,18 @@ func NewPolyKernel(degree int) *PolyKernel {
 
 // InnerProduct computes the inner product through a kernel trick
 // K(x, y) = (x^T y + 1)^d
-func (p *PolyKernel) InnerProduct(vectorX *mat64.Dense, vectorY *mat64.Dense) float64 {
+func (p *PolyKernel) InnerProduct(vectorX *mat.Dense, vectorY *mat.Dense) float64 {
 	subVectorX := vectorX.ColView(0)
 	subVectorY := vectorY.ColView(0)
-	result := mat64.Dot(subVectorX, subVectorY)
+	result := mat.Dot(subVectorX, subVectorY)
 	result = math.Pow(result+1, float64(p.degree))
 
 	return result
 }
 
 // Distance computes distance under the polynomial kernel (maybe not needed?)
-func (p *PolyKernel) Distance(vectorX *mat64.Dense, vectorY *mat64.Dense) float64 {
-	subVector := mat64.NewDense(0, 0, nil)
+func (p *PolyKernel) Distance(vectorX *mat.Dense, vectorY *mat.Dense) float64 {
+	subVector := mat.NewDense(0, 0, nil)
 	subVector.Sub(vectorX, vectorY)
 	result := p.InnerProduct(subVector, subVector)
 
diff --git a/metrics/pairwise/poly_kernel_test.go b/metrics/pairwise/poly_kernel_test.go
index f68e142..8693d79 100644
--- a/metrics/pairwise/poly_kernel_test.go
+++ b/metrics/pairwise/poly_kernel_test.go
@@ -3,17 +3,17 @@ package pairwise
 import (
 	"testing"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	. "github.com/smartystreets/goconvey/convey"
 )
 
 func TestPolyKernel(t *testing.T) {
-	var vectorX, vectorY *mat64.Dense
+	var vectorX, vectorY *mat.Dense
 	polyKernel := NewPolyKernel(3)
 
 	Convey("Given two vectors", t, func() {
-		vectorX = mat64.NewDense(3, 1, []float64{1, 2, 3})
-		vectorY = mat64.NewDense(3, 1, []float64{2, 4, 5})
+		vectorX = mat.NewDense(3, 1, []float64{1, 2, 3})
+		vectorY = mat.NewDense(3, 1, []float64{2, 4, 5})
 
 		Convey("When doing inner product", func() {
 			result := polyKernel.InnerProduct(vectorX, vectorY)
diff --git a/metrics/pairwise/rbf_kernel.go b/metrics/pairwise/rbf_kernel.go
index 1f3fc14..f3d94fa 100644
--- a/metrics/pairwise/rbf_kernel.go
+++ b/metrics/pairwise/rbf_kernel.go
@@ -3,7 +3,7 @@ package pairwise
 import (
 	"math"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type RBFKernel struct {
@@ -17,7 +17,7 @@ func NewRBFKernel(gamma float64) *RBFKernel {
 
 // InnerProduct computes the inner product through a kernel trick
 // K(x, y) = exp(-gamma * ||x - y||^2)
-func (r *RBFKernel) InnerProduct(vectorX *mat64.Dense, vectorY *mat64.Dense) float64 {
+func (r *RBFKernel) InnerProduct(vectorX *mat.Dense, vectorY *mat.Dense) float64 {
 	euclidean := NewEuclidean()
 	distance := euclidean.Distance(vectorX, vectorY)
 
diff --git a/metrics/pairwise/rbf_kernel_test.go b/metrics/pairwise/rbf_kernel_test.go
index 585c5bd..27abd7b 100644
--- a/metrics/pairwise/rbf_kernel_test.go
+++ b/metrics/pairwise/rbf_kernel_test.go
@@ -3,17 +3,17 @@ package pairwise
 import (
 	"testing"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	. "github.com/smartystreets/goconvey/convey"
 )
 
 func TestRBFKernel(t *testing.T) {
-	var vectorX, vectorY *mat64.Dense
+	var vectorX, vectorY *mat.Dense
 	rbfKernel := NewRBFKernel(0.1)
 
 	Convey("Given two vectors", t, func() {
-		vectorX = mat64.NewDense(3, 1, []float64{1, 2, 3})
-		vectorY = mat64.NewDense(3, 1, []float64{2, 4, 5})
+		vectorX = mat.NewDense(3, 1, []float64{1, 2, 3})
+		vectorY = mat.NewDense(3, 1, []float64{2, 4, 5})
 
 		Convey("When doing inner product", func() {
 			result := rbfKernel.InnerProduct(vectorX, vectorY)
diff --git a/naive/bernoulli_nb.go b/naive/bernoulli_nb.go
index 8fb75a3..b128582 100644
--- a/naive/bernoulli_nb.go
+++ b/naive/bernoulli_nb.go
@@ -212,7 +212,7 @@ func (nb *BernoulliNBClassifier) Fit(X base.FixedDataGrid) {
 	docsContainingTerm := make(map[string][]int)
 
 	// This algorithm could be vectorized after binarizing the data
-	// matrix. Since mat64 doesn't have this function, a iterative
+	// matrix. Since mat doesn't have this function, a iterative
 	// version is used.
 	X.MapOverRows(featAttrSpecs, func(docVector [][]byte, r int) (bool, error) {
 		class := base.GetClass(X, r)
diff --git a/neural/layered.go b/neural/layered.go
index fa2c799..3084fb8 100644
--- a/neural/layered.go
+++ b/neural/layered.go
@@ -2,7 +2,7 @@ package neural
 
 import (
 	"fmt"
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	"github.com/sjwhitworth/golearn/base"
 	"github.com/sjwhitworth/golearn/filters"
 	"math"
@@ -105,7 +105,7 @@ func (m *MultiLayerNet) Predict(X base.FixedDataGrid) base.FixedDataGrid {
 	}
 
 	// Create the activation vector
-	a := mat64.NewDense(m.network.size, 1, make([]float64, m.network.size))
+	a := mat.NewDense(m.network.size, 1, make([]float64, m.network.size))
 
 	// Resolve the input AttributeSpecs
 	inputAs := base.ResolveAttributes(insts, inputAttrs)
@@ -277,9 +277,9 @@ func (m *MultiLayerNet) Fit(X base.FixedDataGrid) {
 	}
 
 	// Create the training activation vector
-	trainVec := mat64.NewDense(size, 1, make([]float64, size))
+	trainVec := mat.NewDense(size, 1, make([]float64, size))
 	// Create the error vector
-	errVec := mat64.NewDense(size, 1, make([]float64, size))
+	errVec := mat.NewDense(size, 1, make([]float64, size))
 
 	// Resolve training AttributeSpecs
 	trainAs := base.ResolveAllAttributes(insts)
@@ -322,7 +322,7 @@ func (m *MultiLayerNet) Fit(X base.FixedDataGrid) {
 			}
 
 			// Update total error
-			totalError += math.Abs(mat64.Sum(errVec))
+			totalError += math.Abs(mat.Sum(errVec))
 
 			// Back-propagate the error
 			b := m.network.Error(trainVec, errVec, totalLayers)
diff --git a/neural/layered_test.go b/neural/layered_test.go
index 0a4c496..ac50b14 100644
--- a/neural/layered_test.go
+++ b/neural/layered_test.go
@@ -1,7 +1,7 @@
 package neural
 
 import (
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	"github.com/sjwhitworth/golearn/base"
 	. "github.com/smartystreets/goconvey/convey"
 	"testing"
@@ -129,7 +129,7 @@ func TestLayeredXOR(t *testing.T) {
 				}
 
 			})
-			out := mat64.NewDense(6, 1, []float64{1.0, 0.0, 0.0, 0.0, 0.0, 0.0})
+			out := mat.NewDense(6, 1, []float64{1.0, 0.0, 0.0, 0.0, 0.0, 0.0})
 			net.network.Activate(out, 2)
 			So(out.At(5, 0), ShouldAlmostEqual, 1.0, 0.1)
 
@@ -160,7 +160,7 @@ func TestLayeredXORInline(t *testing.T) {
 
 	Convey("Given an inline XOR dataset...", t, func() {
 
-		data := mat64.NewDense(4, 3, []float64{
+		data := mat.NewDense(4, 3, []float64{
 			1, 0, 1,
 			0, 1, 1,
 			0, 0, 0,
@@ -186,7 +186,7 @@ func TestLayeredXORInline(t *testing.T) {
 				}
 
 			})
-			out := mat64.NewDense(6, 1, []float64{1.0, 0.0, 0.0, 0.0, 0.0, 0.0})
+			out := mat.NewDense(6, 1, []float64{1.0, 0.0, 0.0, 0.0, 0.0, 0.0})
 			net.network.Activate(out, 2)
 			So(out.At(5, 0), ShouldAlmostEqual, 1.0, 0.1)
 
diff --git a/neural/network.go b/neural/network.go
index b5571ed..c00b863 100644
--- a/neural/network.go
+++ b/neural/network.go
@@ -3,7 +3,7 @@ package neural
 import (
 	"bytes"
 	"fmt"
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	"math"
 )
 
@@ -11,8 +11,8 @@ import (
 // Weights are stored in a dense matrix, each can have its own
 // NeuralFunction.
 type Network struct {
-	origWeights *mat64.Dense
-	weights     *mat64.Dense     // n * n
+	origWeights *mat.Dense
+	weights     *mat.Dense     // n * n
 	biases      []float64        // n for each neuron
 	funcs       []NeuralFunction // for each neuron
 	size        int
@@ -27,7 +27,7 @@ type Network struct {
 // connected to themselves for propagation.
 func NewNetwork(size int, input int, f NeuralFunction) *Network {
 	ret := new(Network)
-	ret.weights = mat64.NewDense(size, size, make([]float64, size*size))
+	ret.weights = mat.NewDense(size, size, make([]float64, size*size))
 	ret.biases = make([]float64, size)
 	ret.funcs = make([]NeuralFunction, size)
 	ret.size = size
@@ -104,7 +104,7 @@ func (n *Network) GetBias(node int) float64 {
 // should be set to the number of layers.
 //
 // This function overwrites whatever's stored in its first argument.
-func (n *Network) Activate(with *mat64.Dense, maxIterations int) {
+func (n *Network) Activate(with *mat.Dense, maxIterations int) {
 
 	// Add bias and feed to activation
 	biasFunc := func(r, c int, v float64) float64 {
@@ -114,7 +114,7 @@ func (n *Network) Activate(with *mat64.Dense, maxIterations int) {
 		return n.funcs[r].Forward(v)
 	}
 
-	tmp := new(mat64.Dense)
+	tmp := new(mat.Dense)
 	tmp.Clone(with)
 
 	// Main loop
@@ -128,10 +128,10 @@ func (n *Network) Activate(with *mat64.Dense, maxIterations int) {
 // UpdateWeights takes an output size * 1 output vector and a size * 1
 // back-propagated error vector, as well as a learnRate and updates
 // the internal weights matrix.
-func (n *Network) UpdateWeights(out, err *mat64.Dense, learnRate float64) {
+func (n *Network) UpdateWeights(out, err *mat.Dense, learnRate float64) {
 
 	if n.origWeights == nil {
-		n.origWeights = mat64.DenseCopyOf(n.weights)
+		n.origWeights = mat.DenseCopyOf(n.weights)
 	}
 
 	// Multiply that by the learning rate
@@ -151,7 +151,7 @@ func (n *Network) UpdateWeights(out, err *mat64.Dense, learnRate float64) {
 
 // UpdateBias computes B = B + l.E and updates the bias weights
 // from a size * 1 back-propagated error vector.
-func (n *Network) UpdateBias(err *mat64.Dense, learnRate float64) {
+func (n *Network) UpdateBias(err *mat.Dense, learnRate float64) {
 
 	for i, b := range n.biases {
 		if i < n.input {
@@ -172,11 +172,11 @@ func (n *Network) UpdateBias(err *mat64.Dense, learnRate float64) {
 //
 // If the network is conceptually organised into n layers, maxIterations
 // should be set to n.
-func (n *Network) Error(outArg, errArg *mat64.Dense, maxIterations int) *mat64.Dense {
+func (n *Network) Error(outArg, errArg *mat.Dense, maxIterations int) *mat.Dense {
 
 	// Copy the arguments
-	out := mat64.DenseCopyOf(outArg)
-	err := mat64.DenseCopyOf(errArg)
+	out := mat.DenseCopyOf(outArg)
+	err := mat.DenseCopyOf(errArg)
 
 	// err should be the difference between observed and expected
 	// for observation nodes only (everything else should be zero)
@@ -187,7 +187,7 @@ func (n *Network) Error(outArg, errArg *mat64.Dense, maxIterations int) *mat64.D
 		panic("Unsupported output size")
 	}
 
-	ret := mat64.NewDense(outRows, 1, make([]float64, outRows))
+	ret := mat.NewDense(outRows, 1, make([]float64, outRows))
 
 	// Do differential calculation
 	diffFunc := func(r, c int, v float64) float64 {
@@ -196,7 +196,7 @@ func (n *Network) Error(outArg, errArg *mat64.Dense, maxIterations int) *mat64.D
 	out.Apply(diffFunc, out)
 
 	// Transpose weights matrix
-	reverseWeights := mat64.DenseCopyOf(n.weights)
+	reverseWeights := mat.DenseCopyOf(n.weights)
 	reverseWeights.Clone(n.weights.T())
 
 	// We only need a certain number of passes
diff --git a/neural/network_test.go b/neural/network_test.go
index fee5a29..d59963f 100644
--- a/neural/network_test.go
+++ b/neural/network_test.go
@@ -1,7 +1,7 @@
 package neural
 
 import (
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	. "github.com/smartystreets/goconvey/convey"
 	"testing"
 )
@@ -29,7 +29,7 @@ func TestNetworkWith1Layer(t *testing.T) {
 
 		// Create the Activation vector
 		// NewDense is rows then columns
-		a := mat64.NewDense(6, 1, make([]float64, 6))
+		a := mat.NewDense(6, 1, make([]float64, 6))
 		// Set is rows then columns
 		a.Set(0, 0, 1)
 		a.Set(2, 0, 1)
@@ -42,7 +42,7 @@ func TestNetworkWith1Layer(t *testing.T) {
 			So(a.At(3, 0), ShouldAlmostEqual, 0.332, 0.01)
 
 			// Set the observed error on the output node
-			e := mat64.NewDense(6, 1, make([]float64, 6))
+			e := mat.NewDense(6, 1, make([]float64, 6))
 			e.Set(5, 0, 1.0-a.At(5, 0))
 
 			// Run back-propagated error
diff --git a/neural/neural.go b/neural/neural.go
index a3e9b3b..00301bd 100644
--- a/neural/neural.go
+++ b/neural/neural.go
@@ -2,7 +2,7 @@
 package neural
 
 import (
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type ActivationFunction func(float64) float64
@@ -16,4 +16,4 @@ type NeuralFunction struct {
 
 // LayerFuncs are vectorised layer value transformation functions
 // (e.g. sigmoid). They must operate in-place.
-type LayerFunc func(*mat64.Dense)
+type LayerFunc func(*mat.Dense)
diff --git a/pca/pca.go b/pca/pca.go
index c268157..c6c6837 100644
--- a/pca/pca.go
+++ b/pca/pca.go
@@ -2,13 +2,12 @@
 package pca
 
 import (
-	"github.com/gonum/matrix"
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type PCA struct {
 	Num_components int
-	svd            *mat64.SVD
+	svd            *mat.SVD
 }
 
 // Number of components. 0 - by default, use number of features as number of components
@@ -18,19 +17,19 @@ func NewPCA(num_components int) *PCA {
 
 // Fit PCA model and transform data
 // Need return is base.FixedDataGrid
-func (pca *PCA) FitTransform(X *mat64.Dense) *mat64.Dense {
+func (pca *PCA) FitTransform(X *mat.Dense) *mat.Dense {
 	return pca.Fit(X).Transform(X)
 }
 
 // Fit PCA model
-func (pca *PCA) Fit(X *mat64.Dense) *PCA {
+func (pca *PCA) Fit(X *mat.Dense) *PCA {
 	// Mean to input data
 	M := mean(X)
 	X = matrixSubVector(X, M)
 
 	// Get SVD decomposition from data
-	pca.svd = &mat64.SVD{}
-	ok := pca.svd.Factorize(X, matrix.SVDThin)
+	pca.svd = &mat.SVD{}
+	ok := pca.svd.Factorize(X, mat.SVDThin)
 	if !ok {
 		panic("Unable to factorize")
 	}
@@ -42,41 +41,41 @@ func (pca *PCA) Fit(X *mat64.Dense) *PCA {
 }
 
 // Need return is base.FixedDataGrid
-func (pca *PCA) Transform(X *mat64.Dense) *mat64.Dense {
+func (pca *PCA) Transform(X *mat.Dense) *mat.Dense {
 	if pca.svd == nil {
 		panic("You should to fit PCA model first")
 	}
 
 	num_samples, num_features := X.Dims()
 
-	vTemp := new(mat64.Dense)
-	vTemp.VFromSVD(pca.svd)
+	vTemp := new(mat.Dense)
+	pca.svd.VTo(vTemp)
 	//Compute to full data
 	if pca.Num_components == 0 || pca.Num_components > num_features {
 		return compute(X, vTemp)
 	}
 
 	X = compute(X, vTemp)
-	result := mat64.NewDense(num_samples, pca.Num_components, nil)
-	result.Copy(X.View(0, 0, num_samples, pca.Num_components))
+	result := mat.NewDense(num_samples, pca.Num_components, nil)
+	result.Copy(X)
 	return result
 }
 
 //Helpful private functions
 
 //Compute mean of the columns of input matrix
-func mean(matrix *mat64.Dense) *mat64.Dense {
+func mean(matrix *mat.Dense) *mat.Dense {
 	rows, cols := matrix.Dims()
 	meanVector := make([]float64, cols)
 	for i := 0; i < cols; i++ {
-		sum := mat64.Sum(matrix.ColView(i))
+		sum := mat.Sum(matrix.ColView(i))
 		meanVector[i] = sum / float64(rows)
 	}
-	return mat64.NewDense(1, cols, meanVector)
+	return mat.NewDense(1, cols, meanVector)
 }
 
 // After computing of mean, compute: X(input matrix)  - X(mean vector)
-func matrixSubVector(mat, vec *mat64.Dense) *mat64.Dense {
+func matrixSubVector(mat, vec *mat.Dense) *mat.Dense {
 	rowsm, colsm := mat.Dims()
 	_, colsv := vec.Dims()
 	if colsv != colsm {
@@ -91,8 +90,8 @@ func matrixSubVector(mat, vec *mat64.Dense) *mat64.Dense {
 }
 
 //Multiplication of X(input data) and V(from SVD)
-func compute(X, Y mat64.Matrix) *mat64.Dense {
-	var ret mat64.Dense
+func compute(X, Y mat.Matrix) *mat.Dense {
+	var ret mat.Dense
 	ret.Mul(X, Y)
 	return &ret
 }
diff --git a/pca/pca_test.go b/pca/pca_test.go
index 20595d1..ce0801b 100644
--- a/pca/pca_test.go
+++ b/pca/pca_test.go
@@ -3,13 +3,13 @@ package pca
 import (
 	"testing"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 	. "github.com/smartystreets/goconvey/convey"
 )
 
 func TestPCAWithZeroComponents(t *testing.T) {
 	Convey("Set to pca 0 components with first matrix", t, func() {
-		X1 := mat64.NewDense(3, 7, []float64{6, 5, 4, 3, 8, 2, 9, 5, 1, 10, 2, 3, 8, 7, 5, 14, 2, 3, 6, 3, 2})
+		X1 := mat.NewDense(3, 7, []float64{6, 5, 4, 3, 8, 2, 9, 5, 1, 10, 2, 3, 8, 7, 5, 14, 2, 3, 6, 3, 2})
 		pca := NewPCA(0)
 		rows, cols := pca.FitTransform(X1).Dims()
 		So(rows, ShouldEqual, 3)
@@ -17,7 +17,7 @@ func TestPCAWithZeroComponents(t *testing.T) {
 	})
 
 	Convey("Set to pca 0 components with second matrix", t, func() {
-		X1 := mat64.NewDense(10, 5, []float64{
+		X1 := mat.NewDense(10, 5, []float64{
 			0.52984892, 0.1141001, 0.91599294, 0.9574267, 0.15361222,
 			0.07057588, 0.46371013, 0.73091854, 0.84641034, 0.08122213,
 			0.96221946, 0.60367214, 0.69851546, 0.91965564, 0.27040597,
@@ -37,7 +37,7 @@ func TestPCAWithZeroComponents(t *testing.T) {
 
 func TestPCAWithNComponents(t *testing.T) {
 	Convey("Set to pca 3 components with 5x5 matrix", t, func() {
-		X := mat64.NewDense(5, 5, []float64{
+		X := mat.NewDense(5, 5, []float64{
 			0.23030838, 0.05669317, 0.3187813, 0.34455114, 0.98062806,
 			0.38995469, 0.2996771, 0.99043575, 0.04443827, 0.99527955,
 			0.27266308, 0.14068906, 0.46999473, 0.03296131, 0.90855405,
@@ -50,7 +50,7 @@ func TestPCAWithNComponents(t *testing.T) {
 	})
 
 	Convey("Set to pca 2 components with 3x5 matrix", t, func() {
-		X := mat64.NewDense(3, 5, []float64{
+		X := mat.NewDense(3, 5, []float64{
 			0.12294845, 0.55170713, 0.67572832, 0.60615516, 0.38184551,
 			0.93486821, 0.15120374, 0.89760169, 0.74715672, 0.81373931,
 			0.42821569, 0.47457753, 0.18960954, 0.42466159, 0.34166049})
@@ -63,7 +63,7 @@ func TestPCAWithNComponents(t *testing.T) {
 
 func TestPCAFitAndTransformSeparately(t *testing.T) {
 	Convey("Set to pca 3 components with 5x5 matrix", t, func() {
-		X := mat64.NewDense(5, 5, []float64{
+		X := mat.NewDense(5, 5, []float64{
 			0.23030838, 0.05669317, 0.3187813, 0.34455114, 0.98062806,
 			0.38995469, 0.2996771, 0.99043575, 0.04443827, 0.99527955,
 			0.27266308, 0.14068906, 0.46999473, 0.03296131, 0.90855405,
diff --git a/utilities/utilities.go b/utilities/utilities.go
index f7e5058..d2a9c35 100644
--- a/utilities/utilities.go
+++ b/utilities/utilities.go
@@ -4,7 +4,7 @@ package utilities
 import (
 	"sort"
 
-	"github.com/gonum/matrix/mat64"
+	"gonum.org/v1/gonum/mat"
 )
 
 type sortedIntMap struct {
@@ -37,11 +37,12 @@ func SortIntMap(m map[int]float64) []int {
 	return sm.s
 }
 
-func FloatsToMatrix(floats []float64) *mat64.Dense {
-	return mat64.NewDense(1, len(floats), floats)
+func FloatsToMatrix(floats []float64) *mat.Dense {
+	return mat.NewDense(1, len(floats), floats)
 }
 
-func VectorToMatrix(vector *mat64.Vector) *mat64.Dense {
-	vec := vector.RawVector()
-	return mat64.NewDense(1, len(vec.Data), vec.Data)
+func VectorToMatrix(vector mat.Vector) *mat.Dense {
+	denseCopy := mat.VecDenseCopyOf(vector)
+	vec := denseCopy.RawVector()
+	return mat.NewDense(1, len(vec.Data), vec.Data)
 }