C#. Матрицы: Умножение на скаляр. Инвертирование, транспонирование, детерминант. Сложение, вычитание и умножение

Для реализации множественной линейной регрессии потребовалось поработать с матрицами — а именно, научиться их перемножать, транспонировать и инвертировать. Поскольку ничего такого раньше делать не умел, решил самостоятельно реализовать это в библиотеке машинного обучения, которую сейчас пишу (такой вот способ изучения машинного обучения с нуля придумал для себя).

Библиотека находится на github ВОТ ТУТ. Пока что всё в ней написано на c#. Может быть, еще добавится что-то функциональное — но никак не могу выбрать между F# и всем остальным… По уму-то — нужен F#, так как всё же это у меня .NET. Но вот F#, кажется, сейчас уже никому не нужен, поэтому не знаю, какой язык осваивать.

Итак, создаем матрицу. Конструктор принимает двухмерный массив double — собственно, саму матрицу:

double[,] a = new double[3, 3];
            sqvArr[0, 0] = -1;
            sqvArr[0, 1] = -2;
            sqvArr[0, 2] = 2;
            sqvArr[1, 0] = 2;
            sqvArr[1, 1] = 1;
            sqvArr[1, 2] = 1;
            sqvArr[2, 0] = 3;
            sqvArr[2, 1] = 4;
            sqvArr[2, 2] = 5;

            Matrix matrixA = new Matrix(a);

Также матрицу можно загрузить из файла, указав в методе разделить для столбцов в файле и путь к файлу:

Matrix mFromFile = Matrix.GetMatrixFromTXT("data\\regress_data.txt", '\t');

Если нужно взять часть матрицы (допустим, взять только независимые переменные, или только зависимую (ые), или часть строк), то используется метод, который принимает 2 массива целочисленных значений. 1-ый аргумент — индексы строк, второй — индексы столбцов. В примере ниже берем все строки:

int[] rows = Enumerable.Range(0, mFromFile.matrixBase.GetLength(0))
                        .Select(i => i)
                        .ToArray();
            Matrix bVector = mlr.GetBCoefficientsForMatrix(mFromFile.GetMatrixPart(rows, new int[] { 1, 2, 3 }), 
                mFromFile.GetMatrixPart(rows, new int[] { 0 }));

Для сложения матриц используется статистический метод SumMatrices, принимающий 2 аргумента — суммируемые матрицы и возвращающий их сумму (матрицу, само собой):

Matrix c = Matrix.SumMatrices(matrixA, matrixB);

Вычитаются аналогично:

Matrix c = Matrix.SubstractMatrices(matrixA, matrixB);

Умножение матриц:

Matrix c = Matrix.MultiplyMatrices(matrixA, matrixB);

Умножение матрицы на скаляр уже сделано не статичным:

Matrix c = matrixA.MultiplyToScalar(5);

Транспонирование матрицы:

matrixA = matrixA.Transpose();

Детерминант матрицы:

double determinant = matrixA.GetDeterminant();

Инверсия матрицы:

Matrix invertedSkv = matrixSkv.Invert();

Еще есть переопределенный ToString(). Собственно, если столбцов немного — смотрится нормально.

Самое сложное было реализовать детерминант и инверсию. Вот реализация GetDeterminant():

public double GetDeterminant()
        {
            if (this.matrixBase.GetLength(0) != this.matrixBase.GetLength(1))
                throw new InvalidOperationException("Matrices should be squared");
            
            //if it is 2x2 matrix
            if(this.matrixBase.GetLength(0) == 2)
            {
                return this.matrixBase[0, 0] * this.matrixBase[1, 1]
                    - this.matrixBase[0, 1] * this.matrixBase[1, 0];
            }



            //new matrix with new columns (we need the same number of columns - 1) for the same values 
            //from beginning of matrix
            double[,] extendedMatrixBase = new double[this.matrixBase.GetLength(0), 
                this.matrixBase.GetLength(1) + this.matrixBase.GetLength(1) - 1];

            //filling part of new matrix with the same values as in this object matrix
            for(int i = 0; i < this.matrixBase.GetLength(0); i++)
            {
                for(int j = 0; j < this.matrixBase.GetLength(1); j++)
                {
                    extendedMatrixBase[i, j] = this.matrixBase[i, j];
                }
            }

            //filling new columns with values the same as at beginning of matrix
            for(int i = 0; i < this.matrixBase.GetLength(0); i++)
            {
                for(int j = this.matrixBase.GetLength(1); j < extendedMatrixBase.GetLength(1); j++)
                {
                    extendedMatrixBase[i, j] = this.matrixBase[i, j - this.matrixBase.GetLength(0)];
                }
            }

            Matrix extendedMatrix = new Matrix(extendedMatrixBase);

            //calculating determinant
            double determinant = 0.0;

            //summing
            double sumsResult = 0.0;
            for(int i = 0; i < this.matrixBase.GetLength(0); i++)
            {
                int row = 0;
                int column = i;

                double diagonalResult = 1.0;
                for (int j = 0; j < this.matrixBase.GetLength(0); j++, row++, column++)
                {
                    diagonalResult *= extendedMatrix.matrixBase[row,column];
                }
                sumsResult += diagonalResult;
            }

            //substracting
            double substractionsResult = 0.0;
            for (int i = 0; i < this.matrixBase.GetLength(0); i++)
            {
                int row = 0;
                int column = extendedMatrix.matrixBase.GetLength(1) - 1 - i;

                double diagonalResult = 1.0;
                for(int j = 0; j < this.matrixBase.GetLength(0); j++, row++, column--)
                {
                    diagonalResult *= extendedMatrix.matrixBase[row, column];
                }
                substractionsResult += diagonalResult;
            }

            determinant = sumsResult - substractionsResult;

            return determinant;
        }

И реализация инверсии:

public Matrix Invert()
        {
            Matrix invertedMatrix;

            if (this.isInvertable())
            {
                //Getting matrix of minors
                Matrix matrixOfMinors = new Matrix(new double[matrixBase.GetLength(0), matrixBase.GetLength(1)]);

                for (int row = 0; row < matrixBase.GetLength(0); row++)
                {
                    for(int column= 0; column < matrixBase.GetLength(1); column++)
                    {
                        double[,] subMArrBase = new double[matrixBase.GetLength(0)-1,matrixBase.GetLength(1)-1];
                        Matrix subMatrix = new Matrix(subMArrBase);

                        int subMRow = 0;
                        for(int r = 0; r < matrixBase.GetLength(0); r++)
                        {
                            if (row == r)
                                continue;
                            int subMColumn = 0;
                            for (int c = 0; c < matrixBase.GetLength(1); c++)
                            {
                                if (column == c)
                                    continue;

                                subMArrBase[subMRow, subMColumn] = this.matrixBase[r, c];
                                subMColumn++;
                            }
                            subMRow++;
                        }

                        matrixOfMinors.matrixBase[row, column] = subMatrix.GetDeterminant();
                    }
                }

                //creating cofactor matrix
                int rowStartMultiplier;
                for(int i = 0; i < matrixOfMinors.matrixBase.GetLength(0); i++)
                {
                    if (i % 2 == 0)
                        rowStartMultiplier = 1;
                    else
                        rowStartMultiplier = -1;

                    for(int j = 0; j < matrixOfMinors.matrixBase.GetLength(1); j++)
                    {
                        matrixOfMinors.matrixBase[i, j] *= rowStartMultiplier;
                        //swap
                        rowStartMultiplier *= -1;
                    }
                }

                //determinant of THIS matrix
                double determinant = this.GetDeterminant();
                Matrix transposed = matrixOfMinors.Transpose();
                invertedMatrix = transposed.MultiplyToScalar(1 / determinant);

            }
            else
                throw new InvalidOperationException("Matrix is not invertable");

            return invertedMatrix;
        }

В следующем посте покажу, как написал множественную линейную регрессию.

 

Добавить комментарий