using System;
using System.Collections.Generic;
using System.ComponentModel;
using System.Diagnostics;
using System.Drawing;
using System.Linq;
using System.Runtime.InteropServices;
using System.Security.Cryptography;
using System.Text;
using System.Threading;
using System.Threading.Tasks;
using System.Transactions;
namespace MulLinearRegression
{
public enum NormType
{
L0,
L1,
L2,
LP,
P_INF,
M_INF,
INF,
F
}
public class Matrix
{
public static double[,] Transpose(double[,] iMatrix)
{
int row = iMatrix.GetLength(0);
int column = iMatrix.GetLength(1);
double[,] TempMatrix = new double[row, column];
double[,] iMatrixT = new double[column, row];
for (int i = 0; i < row; i++)
{
for (int j = 0; j < column; j++)
{
TempMatrix[i, j] = iMatrix[i, j];
}
}
for (int i = 0; i < column; i++)
{
for (int j = 0; j < row; j++)
{
iMatrixT[i, j] = TempMatrix[j, i];
}
}
return iMatrixT;
}
public static double[,] Athwart(double[,] iMatrix)
{
int i = 0;
int row = iMatrix.GetLength(0);
double[,] MatrixZwei = new double[row, row * 2];
double[,] iMatrixInv = new double[row, row];
for (i = 0; i < row; i++)
{
for (int j = 0; j < row; j++)
{
MatrixZwei[i, j] = iMatrix[i, j];
}
}
for (i = 0; i < row; i++)
{
for (int j = row; j < row * 2; j++)
{
MatrixZwei[i, j] = 0;
if (i + row == j)
MatrixZwei[i, j] = 1;
}
}
for (i = 0; i < row; i++)
{
if (MatrixZwei[i, i] != 0)
{
double intTemp = MatrixZwei[i, i];
for (int j = 0; j < row * 2; j++)
{
MatrixZwei[i, j] = MatrixZwei[i, j] / intTemp;
}
}
for (int j = 0; j < row; j++)
{
if (j == i)
continue;
double intTemp = MatrixZwei[j, i];
for (int k = 0; k < row * 2; k++)
{
MatrixZwei[j, k] = MatrixZwei[j, k] - MatrixZwei[i, k] * intTemp;
}
}
}
for (i = 0; i < row; i++)
{
for (int j = 0; j < row; j++)
{
iMatrixInv[i, j] = MatrixZwei[i, j + row];
}
}
return iMatrixInv;
}
public static double[,] AddMatrix(double[,] MatrixEin, double[,] MatrixZwei)
{
double[,] MatrixResult = new double[MatrixEin.GetLength(0), MatrixZwei.GetLength(1)];
for (int i = 0; i < MatrixEin.GetLength(0); i++)
for (int j = 0; j < MatrixZwei.GetLength(1); j++)
MatrixResult[i, j] = MatrixEin[i, j] + MatrixZwei[i, j];
return MatrixResult;
}
public static double[,] MultiplyMatrix(double[,] MatrixEin, double[,] MatrixZwei)
{
double[,] MatrixResult = new double[MatrixEin.GetLength(0), MatrixZwei.GetLength(1)];
for (int i = 0; i < MatrixEin.GetLength(0); i++)
{
for (int j = 0; j < MatrixZwei.GetLength(1); j++)
{
for (int k = 0; k < MatrixEin.GetLength(1); k++)
{
MatrixResult[i, j] += MatrixEin[i, k] * MatrixZwei[k, j];
}
}
}
return MatrixResult;
}
public static double[,] HadamardProduct(double[,] MatrixEin, double[,] MatrixZwei)
{
int Row = MatrixEin.GetLength(0);
int Col=MatrixEin.GetLength(1);
double[,] outMatrix = new double[Row, Col];
for(int i = 0; i < Row; i++)
{
for(int j = 0; j < Col; j++)
{
outMatrix[i, j] = MatrixEin[i, j] * MatrixZwei[i, j];
}
}
return outMatrix;
}
public static double ResultDeterminant(double[,] MatrixEin)
{
return MatrixEin[0, 0] * MatrixEin[1, 1] * MatrixEin[2, 2] + MatrixEin[0, 1] * MatrixEin[1, 2] * MatrixEin[2, 0] + MatrixEin[0, 2] * MatrixEin[1, 0] * MatrixEin[2, 1]
- MatrixEin[0, 2] * MatrixEin[1, 1] * MatrixEin[2, 0] - MatrixEin[0, 1] * MatrixEin[1, 0] * MatrixEin[2, 2] - MatrixEin[0, 0] * MatrixEin[1, 2] * MatrixEin[2, 1];
}
public static double[,] UnitMatrix(int Row)
{
double[,] outMatrix = new double[Row, Row];
for (int i = 0; i < Row; i++)
outMatrix[i, i] = 1;
return outMatrix;
}
public static double[,] VectorToDiagonalMatrix(double[,] MatrixEin)
{
int Row = MatrixEin.GetLength(0);
int Col = MatrixEin.GetLength(1);
double[] data = SubDatas(MatrixEin, 0, RowGet: Col > Row);
Row = Row > Col ? Row : Col;
double[,] outMatrix = new double[Row, Row];
for(int i = 0; i < Row; i++)
{
outMatrix[i, i] = data[i];
}
return outMatrix;
}
public static double[,] MatrixCombine(double[,] MatrixEin, double[,] MatrixZwei, bool RowCombine = false)
{
int rowE = MatrixEin.GetLength(0);
int colE = MatrixEin.GetLength(1);
int rowZ = MatrixZwei.GetLength(0);
int colZ = MatrixZwei.GetLength(1);
int dstRow = RowCombine ? rowE + rowZ : rowE;
int dstCol = RowCombine ? colE : colE + colZ;
double[,] outMatrix = new double[dstRow, dstCol];
if (rowE * colE == 0)
{
return (double[,])MatrixZwei.Clone();
}
if (rowZ * colZ == 0)
{
return (double[,])MatrixEin.Clone();
}
if (RowCombine)
{
Array.Copy(MatrixEin, 0, outMatrix, 0, rowE * colE);
Array.Copy(MatrixZwei, 0, outMatrix, rowE * colE, rowZ * colZ);
}
else
{
outMatrix = new double[dstCol, dstRow];
MatrixEin = Transpose(MatrixEin);
MatrixZwei = Transpose(MatrixZwei);
Array.Copy(MatrixEin, 0, outMatrix, 0, rowE * colE);
Array.Copy(MatrixZwei, 0, outMatrix, rowE * colE, rowZ * colZ);
outMatrix = Transpose(outMatrix);
}
return outMatrix;
}
public static double[,] SubMatrix(double[,] MatrixEin,int Start,int len=1,bool RowGet = false)
{
double[,] temp = RowGet? MatrixEin:Transpose(MatrixEin);
int Row=temp.GetLength(0);
int Col=temp.GetLength(1);
double[,] outMatrix = new double[len, Col];
Buffer.BlockCopy(temp,sizeof(double)*Start*Col, outMatrix, 0, sizeof(double) * Col * len);
return RowGet? outMatrix:Transpose(outMatrix);
}
public static double[,] SubMatrix(double[,] MatrixEin,Point LeftTop,int Width,int Height)
{
double[,] outMatrix;
int uRow = LeftTop.Y;
int lCol = LeftTop.X;
int dRow = uRow + Height;
int rCol = lCol + Width;
outMatrix = SubMatrix(MatrixEin, uRow, Height, RowGet: true);
outMatrix = SubMatrix(outMatrix, lCol, Width, RowGet: false);
return outMatrix;
}
public static double[] SubDatas(double[,] MatrixEin,int Start,int len=1,bool RowGet = false)
{
double[,] temp = RowGet ? (double[,])MatrixEin.Clone() : Transpose(MatrixEin);
int Row = temp.GetLength(0);
int Col = temp.GetLength(1);
double[] data = new double[Col*len];
Buffer.BlockCopy(temp, Start * Col * sizeof(double), data, 0, Col * len * sizeof(double));
return data;
}
public static double[,] VectorGenerate(double[] Data, bool RowVec = false)
{
double[,] outVector = new double[1, Data.Length];
Buffer.BlockCopy(Data, 0, outVector, 0, sizeof(double) * Data.Length);
outVector = RowVec ? outVector : Transpose(outVector);
return outVector;
}
public static double[,] VectorGenerate(int len, double value = 1.0, bool RowVec = false)
{
double[,] array = new double[1, len];
for (int i = 0; i < len; i++)
{
array[0, i] = value;
}
if (RowVec)
{
return array;
}
else
{
return Transpose(array);
}
}
public static double[,] MatrixRemoveVector(double[,] MatrixEin, int Id, bool RowRmv = false)
{
double[,] uMatrix = MultiplyConst(MatrixEin, 1);
uMatrix = RowRmv ? uMatrix : Transpose(uMatrix);
int eRow = uMatrix.GetLength(0);
int eCol = uMatrix.GetLength(1);
int dRow = eRow - 1;
double[,] outMatrix = new double[dRow, eCol];
Array.Copy(uMatrix, 0, outMatrix, 0, Id * eCol);
if (Id + 1 < eRow)
{
Array.Copy(uMatrix, (Id + 1) * eCol, outMatrix, Id * eCol,
(eRow - Id - 1) * eCol);
}
outMatrix = RowRmv ? outMatrix : Transpose(outMatrix);
return outMatrix;
}
public static double Norm(double[,] MatrixEin, NormType type, double P = 0)
{
double norm = double.MinValue;
int Row = MatrixEin.GetLength(0);
int Col = MatrixEin.GetLength(1);
if (Row * Col == Row || Row * Col == Col)
{
double[] data = SubDatas(MatrixEin, 0, RowGet: Col > Row);
switch (type)
{
case NormType.L0:
{
norm = data.ToList().Count(x => x != 0);
break;
}
case NormType.L1:
{
norm = data.Sum(x => Math.Abs(x));
break;
}
case NormType.L2:
{
norm = Math.Sqrt(data.Sum(x => x * x));
break;
}
case NormType.LP:
{
norm = Math.Pow(data.Sum(x => Math.Pow(Math.Abs(x), P)), 1.0 / (double)P);
break;
}
case NormType.P_INF:
{
norm = data.Max(x => Math.Abs(x));
break;
}
case NormType.M_INF:
{
norm = data.Min(x => Math.Abs(x));
break;
}
}
}
else
{
switch (type)
{
case NormType.L1:
{
for (int i = 0; i < Col; i++)
{
norm = Math.Max(norm, SubDatas(MatrixEin, i, RowGet: false).Sum(x => Math.Abs(x)));
}
break;
}
case NormType.INF:
{
for (int i = 0; i < Row; i++)
{
norm = Math.Max(norm, SubDatas(MatrixEin, i, RowGet: true).Sum(x => Math.Abs(x)));
}
break;
}
case NormType.F:
{
norm = 0;
for (int i = 0; i < Row; i++)
{
norm += SubDatas(MatrixEin, i, RowGet: true).Sum(x => x * x);
}
norm = Math.Sqrt(norm);
break;
}
case NormType.L2:
{
double[,] nnd = Transpose(MatrixEin);
nnd = MultiplyMatrix(nnd, MatrixEin);
double lambda = MatrixEigenvalue(nnd, ref nnd, Type: 0);
norm = Math.Sqrt(lambda);
break;
}
}
}
return norm;
}
public static int MatrixTypeJudge(double[,] Datas)
{
int Row = Datas.GetLength(0);
int Col = Datas.GetLength(1);
if (Row == 1 && Col == 1)
{
return 0;
}
else if (Row * Col == Row || Row * Col == Col)
{
return 1;
}
else
{
return 2;
}
}
public static void MinMax(double[,] MatrixEin, out double Min, out double Max)
{
Max = double.MinValue;
Min = double.MaxValue;
int Row = MatrixEin.GetLength(0);
int Col = MatrixEin.GetLength(1);
for (int i = 0; i < Row; i++)
{
Max = Math.Max(Max, SubDatas(MatrixEin, i, RowGet: true).Max());
Min = Math.Min(Min, SubDatas(MatrixEin, i, RowGet: true).Min());
}
}
public static double[,] Normalization(double[,] MatrixEin,int type = 0)
{
int Row=MatrixEin.GetLength(0);
int Col = MatrixEin.GetLength(1);
double[,] outMatrix = new double[Row, Col];
switch (type)
{
case 0:
{
double squSum = Math.Pow(Norm(MatrixEin, NormType.F), 2);
for (int i = 0; i < Row; i++)
{
for (int j = 0; j < Col; j++)
{
outMatrix[i, j] = MatrixEin[i, j] * MatrixEin[i, j] / squSum;
}
}
break;
}
case 1:
{
for (int c = 0; c < Col; c++)
{
double colSquSum = Math.Pow(Norm(SubMatrix(MatrixEin, c), NormType.L2), 2);
for (int r = 0; r < Row; r++)
{
outMatrix[r, c] = MatrixEin[r, c] * MatrixEin[r, c] / colSquSum;
}
}
break;
}
case 2:
{
double[] ljj = new double[Col];
double[] avgCol = new double[Col];
for (int c = 0; c < Col; c++)
{
avgCol[c] = SubDatas(MatrixEin, c).Average();
for (int r = 0; r < Row; r++)
{
ljj[c] += Math.Pow(MatrixEin[r, c] - avgCol[c], 2);
}
}
for(int c = 0; c < Col; c++)
{
for(int r = 0; r < Row; r++)
{
outMatrix[r, c] = (MatrixEin[r, c] - avgCol[c]) / Math.Sqrt(ljj[c]);
}
}
break;
}
case 3:
{
for(int c = 0; c < Col; c++)
{
double colAbsMax = Norm(SubMatrix(MatrixEin, c), NormType.P_INF);
for(int r = 0; r < Row; r++)
{
outMatrix[r, c] = MatrixEin[r, c] / colAbsMax;
}
}
break;
}
}
return outMatrix;
}
public static double MatrixEigenvalue(double[,] MatrixEin, ref double[,] Eigenvector, double epsilon = 1e-6, int Type = 0)
{
MatrixEin = Type == 0 ? MatrixEin : Athwart(MatrixEin);
double Eigenvalue = double.MaxValue;
int Row = MatrixEin.GetLength(0);
double[,] cVec = VectorGenerate(Row, 1, RowVec: false);
double[,] uVec = VectorGenerate(Row, 1, RowVec: false);
double diff = double.MaxValue;
while (diff > epsilon)
{
cVec = MultiplyMatrix(MatrixEin, uVec);
uVec = Normalization(cVec, 3);
MinMax(cVec, out double min, out double max);
diff = Math.Abs(max - Eigenvalue);
Eigenvalue = max;
Eigenvector = uVec;
}
return Eigenvalue;
}
public static double[] SymmetricMatricEigenvalud(double[,] MatrixEin, ref double[,] EigenVector, double epsilon = 1e-7, int iterations = 1000)
{
double[,] A = (double[,])MatrixEin.Clone();
int Row = MatrixEin.GetLength(0);
int Col = MatrixEin.GetLength(1);
double[,] Q = UnitMatrix(Row);
double[,] R = UnitMatrix(Row);
double max = double.MaxValue;
int mx = -1, my = -1;
int times = 0;
while (max > epsilon && times < iterations)
{
max = double.MinValue;
for (int i = 0; i < Row; i++)
{
for (int j = i + 1; j < Col; j++)
{
if (Math.Abs(A[i, j]) > max)
{
max = Math.Abs(A[i, j]);
mx = i;
my = j;
}
}
}
times += 1;
double d = (A[mx, mx] - A[my, my]) / (2.0 * A[mx, my]);
double t = d > 0 ? -d + Math.Sqrt(d * d + 1) : d == 0 ? -1 : -d - Math.Sqrt(d * d + 1);
double c = Math.Pow(1 + t * t, -0.5);
double s = t * c;
R = UnitMatrix(Row);
R[mx, mx] = R[my, my] = c;
R[mx, my] = s;
R[my, mx] = -s;
double[,] RT = Transpose(R);
double[,] B = MultiplyMatrix(R, A);
A = MultiplyMatrix(B, RT);
Q = MultiplyMatrix(Q, RT);
}
EigenVector = Q;
double[] Eigenvalue = new double[Col];
for (int i = 0; i < Row; i++)
{
Eigenvalue[i] = A[i, i];
}
return Eigenvalue;
}
public static void MatrixShuffle(double[,] MatrixEin, out double[,] shuffled)
{
int Row = MatrixEin.GetLength(0);
int Col = MatrixEin.GetLength(1);
shuffled = new double[Row, Col];
int[] randIds = Enumerable.Range(0, Row).OrderBy(x => Guid.NewGuid()).Take(Row).ToArray();
for (int i = 0; i < randIds.Length; i++)
{
Buffer.BlockCopy(MatrixEin, randIds[i] * Col * sizeof(double), shuffled, i * Col * sizeof(double), Col * sizeof(double));
}
}
public static double[,] MultiplyConst(double[,] MatrixEin, double Coefficient)
{
int Row = MatrixEin.GetLength(0);
int Col = MatrixEin.GetLength(1);
double[,] outMatrix = new double[Row, Col];
for (int i = 0; i < outMatrix.GetLength(0); i++)
{
for (int j = 0; j < outMatrix.GetLength(1); j++)
{
outMatrix[i, j] = MatrixEin[i, j] * Coefficient;
}
}
return outMatrix;
}
public static double[,] MatrixElementSqrt(double[,] MatrixEin)
{
int Row = MatrixEin.GetLength(0);
int Col = MatrixEin.GetLength(1);
double[,] outMatrix = new double[Row, Col];
for (int i = 0; i < Row; i++)
{
for (int j = 0; j < Col; j++)
{
outMatrix[i, j] = Math.Sqrt(MatrixEin[i, j]);
}
}
return outMatrix;
}
public static double[,] MatrixElementReciprocal(double[,] MatrixEin)
{
int Row = MatrixEin.GetLength(0);
int Col = MatrixEin.GetLength(1);
double[,] outMatrix = new double[Row, Col];
for (int i = 0; i < Row; i++)
{
for (int j = 0; j < Col; j++)
{
outMatrix[i, j] = 1.0 / MatrixEin[i, j];
}
}
return outMatrix;
}
public static double[,] VectorStretch(double[,] Vector,int Size = 1)
{
int Row=Vector.GetLength(0);
int Col=Vector.GetLength(1);
double[,] temp = Row > Col ? Transpose(Vector) : (double[,])Vector.Clone();
double[,] outMatrix = new double[Size, Row > Col ? Row : Col];
for(int i = 0; i < Size; i++)
{
Buffer.BlockCopy(temp, 0, outMatrix, i * sizeof(double) * (Row > Col ? Row : Col), sizeof(double) * (Row > Col ? Row : Col));
}
return Row > Col ? Transpose(outMatrix) : outMatrix;
}
public static double[,] MatrixReshape(double[,] MatrixEin,int Row)
{
int Width = MatrixEin.GetLength(0) * MatrixEin.GetLength(1) / Row;
double[,] outMatrix = new double[Row,Width];
for(int i = 0; i < Row; i++)
{
Buffer.BlockCopy(MatrixEin, i * Width * sizeof(double), outMatrix, i * Width * sizeof(double), sizeof(double) * Width);
}
return outMatrix;
}
public static double[,] MatrixSigmoid(double[,] MatrixEin,bool Grad=false)
{
double[,] outMatrix = new double[MatrixEin.GetLength(0), MatrixEin.GetLength(1)];
if (!Grad)
{
for (int i = 0; i < outMatrix.GetLength(0); i++)
for (int j = 0; j < outMatrix.GetLength(1); j++)
outMatrix[i, j] = 1.0 / (1 + Math.Exp(-1 * MatrixEin[i, j]));
}
else
{
for (int i = 0; i < outMatrix.GetLength(0); i++)
for (int j = 0; j < outMatrix.GetLength(1); j++)
outMatrix[i, j] = MatrixEin[i, j] * (1 - MatrixEin[i, j]);
}
return outMatrix;
}
public static double[,] MatrixTanh(double[,] MatrixEin, bool Grad = false)
{
double[,] outMatrix = new double[MatrixEin.GetLength(0), MatrixEin.GetLength(1)];
if (!Grad)
{
for (int i = 0; i < outMatrix.GetLength(0); i++)
for (int j = 0; j < outMatrix.GetLength(1); j++)
outMatrix[i, j] = (1.0 - Math.Exp(-2 * MatrixEin[i, j])) / (1.0 + Math.Exp(-2 * MatrixEin[i, j]));
}
else
{
for (int i = 0; i < outMatrix.GetLength(0); i++)
for (int j = 0; j < outMatrix.GetLength(1); j++)
outMatrix[i, j] = 1-(MatrixEin[i, j]* MatrixEin[i, j]);
}
return outMatrix;
}
public static double[,] MatrixReLU(double[,] MatrixEin, bool Grad = false)
{
double[,] outMatrix = new double[MatrixEin.GetLength(0), MatrixEin.GetLength(1)];
if (!Grad)
{
for (int i = 0; i < outMatrix.GetLength(0); i++)
for (int j = 0; j < outMatrix.GetLength(1); j++)
outMatrix[i, j] = Math.Max(0, MatrixEin[i, j]);
}
else
{
for (int i = 0; i < outMatrix.GetLength(0); i++)
for (int j = 0; j < outMatrix.GetLength(1); j++)
outMatrix[i, j] = MatrixEin[i, j] > 0 ? 1 : 0;
}
return outMatrix;
}
public static double[,] MatrixLog(double[,] MatrixEin,double basic=2)
{
int row = MatrixEin.GetLength(0);
int col = MatrixEin.GetLength(1);
double[,] outMatrix = new double[row, col];
for(int i = 0; i < row; i++)
{
for(int j = 0; j < col; j++)
{
outMatrix[i, j] =
basic == 2 ? Math.Log2(MatrixEin[i, j]) :
basic == 10 ? Math.Log10(MatrixEin[i, j]) :
Math.Log(MatrixEin[i, j]);
}
}
return outMatrix;
}
}
}