using System;
using System.Collections.Generic;
using System.Linq;
using System.Data;
namespace MulLinearRegression
{
class Regression
{
//α=0.01
static public double[,] F_Dist_1E_2 = {
{
4052 , 4999 , 5403 , 5625 , 5764 , 5859 , 5981 , 6106 , 6234 , 6366 } ,
{
98.49 , 99.01 , 99.17 , 99.25 , 99.3 , 99.33 , 99.36 , 99.42 , 99.46 , 99.5 } ,
{
34.12 , 30.81 , 29.46 , 28.71 , 28.24 , 27.91 , 27.49 , 27.05 , 26.6 , 26.12 } ,
{
21.2 , 18 , 16.69 , 15.98 , 15.52 , 15.21 , 14.8 , 14.37 , 13.93 , 13.46 } ,
{
16.26 , 13.27 , 12.06 , 11.39 , 10.97 , 10.67 , 10.29 , 9.89 , 9.47 , 9.02 } ,
{
13.74 , 10.92 , 9.78 , 9.15 , 8.75 , 8.47 , 8.1 , 7.72 , 7.31 , 6.88 } ,
{
12.25 , 9.55 , 8.45 , 7.85 , 7.46 , 7.19 , 6.84 , 6.47 , 6.07 , 5.65 } ,
{
11.26 , 8.65 , 7.59 , 7.01 , 6.63 , 6.37 , 6.03 , 5.67 , 5.28 , 4.86 } ,
{
10.56 , 8.02 , 6.99 , 6.42 , 6.06 , 5.8 , 5.47 , 5.11 , 4.73 , 4.31 } ,
{
10.04 , 7.56 , 6.55 , 5.99 , 5.64 , 5.39 , 5.06 , 4.71 , 4.33 , 3.91 } ,
{
9.65 , 7.2 , 6.22 , 5.67 , 5.32 , 5.07 , 4.74 , 4.4 , 4.02 , 3.6 } ,
{
9.33 , 6.93 , 5.95 , 5.41 , 5.06 , 4.82 , 4.5 , 4.16 , 3.78 , 3.36 } ,
{
9.07 , 6.7 , 5.74 , 5.2 , 4.86 , 4.62 , 4.3 , 3.96 , 3.59 , 3.16 } ,
{
8.86 , 6.51 , 5.56 , 5.03 , 4.69 , 4.46 , 4.14 , 3.8 , 3.43 , 3 } ,
{
8.68 , 6.36 , 5.42 , 4.89 , 4.56 , 4.32 , 4 , 3.67 , 3.29 , 2.87 } ,
{
8.53 , 6.23 , 5.29 , 4.77 , 4.44 , 4.2 , 3.89 , 3.55 , 3.18 , 2.75 } ,
{
8.4 , 6.11 , 5.18 , 4.67 , 4.34 , 4.1 , 3.79 , 3.45 , 3.08 , 2.65 } ,
{
8.28 , 6.01 , 5.09 , 4.58 , 4.25 , 4.01 , 3.71 , 3.37 , 3 , 2.57 } ,
{
8.18 , 5.93 , 5.01 , 4.5 , 4.17 , 3.94 , 3.63 , 3.3 , 2.92 , 2.49 } ,
{
8.1 , 5.85 , 4.94 , 4.43 , 4.1 , 3.87 , 3.56 , 3.23 , 2.86 , 2.42 } ,
{
8.02 , 5.78 , 4.87 , 4.37 , 4.04 , 3.81 , 3.51 , 3.17 , 2.8 , 2.36 } ,
{
7.94 , 5.72 , 4.82 , 4.31 , 3.99 , 3.76 , 3.45 , 3.12 , 2.75 , 2.31 } ,
{
7.88 , 5.66 , 4.76 , 4.26 , 3.94 , 3.71 , 3.41 , 3.07 , 2.7 , 2.26 } ,
{
7.82 , 5.61 , 4.72 , 4.22 , 3.9 , 3.67 , 3.36 , 3.03 , 2.66 , 2.21 } ,
{
7.77 , 5.57 , 4.68 , 4.18 , 3.86 , 3.63 , 3.32 , 2.99 , 2.62 , 2.17 } ,
{
7.72 , 5.53 , 4.64 , 4.14 , 3.82 , 3.59 , 3.29 , 2.96 , 2.58 , 2.13 } ,
{
7.68 , 5.49 , 4.6 , 4.11 , 3.78 , 3.56 , 3.26 , 2.93 , 2.55 , 2.1 } ,
{
7.64 , 5.45 , 4.57 , 4.07 , 3.75 , 3.53 , 3.23 , 2.9 , 2.52 , 2.06 } ,
{
7.6 , 5.42 , 4.54 , 4.04 , 3.73 , 3.5 , 3.2 , 2.87 , 2.49 , 2.03 } ,
{
7.56 , 5.39 , 4.51 , 4.02 , 3.7 , 3.47 , 3.17 , 2.84 , 2.47 , 2.01 } ,
{
7.31 , 5.18 , 4.31 , 3.83 , 3.51 , 3.29 , 2.99 , 2.66 , 2.29 , 1.8 } ,
{
7.08 , 4.98 , 4.13 , 3.65 , 3.34 , 3.12 , 2.82 , 2.5 , 2.12 , 1.6 } ,
{
6.85 , 4.79 , 3.95 , 3.48 , 3.17 , 2.96 , 2.66 , 2.34 , 1.95 , 1.38 } ,
{
6.64 , 4.6 , 3.78 , 3.32 , 3.02 , 2.8 , 2.51 , 2.18 , 1.79 , 1 }
};
//α=0.025
static public double[,] F_Dist_2p5E_2 = {
{
647.8 , 799.5 , 864.2 , 899.6 , 921.8 , 937.1 , 956.7 , 976.7 , 997.2 , 1018 } ,
{
38.51 , 39 , 39.17 , 39.25 , 39.3 , 39.33 , 39.37 , 39.41 , 39.46 , 39.5 } ,
{
17.44 , 16.04 , 15.44 , 15.1 , 14.88 , 14.73 , 14.54 , 14.34 , 14.12 , 13.9 } ,
{
12.22 , 10.65 , 9.98 , 9.6 , 9.36 , 9.2 , 8.98 , 8.75 , 8.51 , 8.26 } ,
{
10.01 , 8.43 , 7.76 , 7.39 , 7.15 , 6.98 , 6.76 , 6.52 , 6.28 , 6.02 } ,
{
8.81 , 7.26 , 6.6 , 6.23 , 5.99 , 5.82 , 5.6 , 5.37 , 5.12 , 4.85 } ,
{
8.07 , 6.54 , 5.89 , 5.52 , 5.29 , 5.12 , 4.9 , 4.67 , 4.42 , 4.14 } ,
{
7.57 , 6.06 , 5.42 , 5.05 , 4.82 , 4.65 , 4.43 , 4.2 , 3.95 , 3.67 } ,
{
7.21 , 5.71 , 5.08 , 4.72 , 4.48 , 4.32 , 4.1 , 3.87 , 3.61 , 3.33 } ,
{
6.94 , 5.46 , 4.83 , 4.47 , 4.24 , 4.07 , 3.85 , 3.62 , 3.37 , 3.08 } ,
{
6.72 , 5.26 , 4.63 , 4.28 , 4.04 , 3.88 , 3.66 , 3.43 , 3.17 , 2.88 } ,
{
6.55 , 5.1 , 4.47 , 4.12 , 3.89 , 3.73 , 3.51 , 3.28 , 3.02 , 2.72 } ,
{
6.41 , 4.97 , 4.35 , 4 , 3.77 , 3.6 , 3.39 , 3.15 , 2.89 , 2.6 } ,
{
6.3 , 4.86 , 4.24 , 3.89 , 3.66 , 3.5 , 3.29 , 3.05 , 2.79 , 2.49 } ,
{
6.2 , 4.77 , 4.15 , 3.8 , 3.58 , 3.41 , 3.2 , 2.96 , 2.7 , 2.4 } ,
{
6.12 , 4.69 , 4.08 , 3.73 , 3.5 , 3.34 , 3.12 , 2.89 , 2.63 , 2.32 } ,
{
6.04 , 4.62 , 4.01 , 3.66 , 3.44 , 3.28 , 3.06 , 2.82 , 2.56 , 2.25 } ,
{
5.98 , 4.56 , 3.95 , 3.61 , 3.38 , 3.22 , 3.01 , 2.77 , 2.5 , 2.19 } ,
{
5.92 , 4.51 , 3.9 , 3.56 , 3.33 , 3.17 , 2.96 , 2.72 , 2.45 , 2.13 } ,
{
5.87 , 4.46 , 3.86 , 3.51 , 3.29 , 3.13 , 2.91 , 2.68 , 2.41 , 2.09 } ,
{
5.83 , 4.42 , 3.82 , 3.48 , 3.25 , 3.09 , 2.87 , 2.64 , 2.37 , 2.04 } ,
{
5.79 , 4.38 , 3.78 , 3.44 , 3.22 , 3.05 , 2.84 , 2.6 , 2.33 , 2 } ,
{
5.75 , 4.35 , 3.75 , 3.41 , 3.18 , 3.02 , 2.81 , 2.57 , 2.3 , 1.97 } ,
{
5.72 , 4.32 , 3.72 , 3.38 , 3.15 , 2.99 , 2.78 , 2.54 , 2.27 , 1.94 } ,
{
5.69 , 4.29 , 3.69 , 3.35 , 3.13 , 2.97 , 2.75 , 2.51 , 2.24 , 1.91 } ,
{
5.66 , 4.27 , 3.67 , 3.33 , 3.1 , 2.94 , 2.73 , 2.49 , 2.22 , 1.88 } ,
{
5.63 , 4.24 , 3.65 , 3.31 , 3.08 , 2.92 , 2.71 , 2.47 , 2.19 , 1.85 } ,
{
5.61 , 4.22 , 3.63 , 3.29 , 3.06 , 2.9 , 2.69 , 2.45 , 2.17 , 1.83 } ,
{
5.59 , 4.2 , 3.61 , 3.27 , 3.04 , 2.88 , 2.67 , 2.43 , 2.15 , 1.81 } ,
{
5.57 , 4.18 , 3.59 , 3.25 , 3.03 , 2.87 , 2.65 , 2.41 , 2.14 , 1.79 } ,
{
5.42 , 4.05 , 3.46 , 3.13 , 2.9 , 2.74 , 2.53 , 2.29 , 2.01 , 1.64 } ,
{
5.29 , 3.93 , 3.34 , 3.01 , 2.79 , 2.63 , 2.41 , 2.17 , 1.88 , 1.48 } ,
{
5.15 , 3.8 , 3.23 , 2.89 , 2.67 , 2.52 , 2.3 , 2.05 , 1.76 , 1.31 } ,
{
5.02 , 3.69 , 3.12 , 2.79 , 2.57 , 2.41 , 2.19 , 1.94 , 1.64 , 1 }
};
//α=0.05
static public double[,] F_Dist_5E_2 = {
{
161.4 , 199.5 , 215.7 , 224.6 , 230.2 , 234 , 238.9 , 243.9 , 249 , 254.3 } ,
{
18.51 , 19 , 19.16 , 19.25 , 19.3 , 19.33 , 19.37 , 19.41 , 19.45 , 19.5 } ,
{
10.13 , 9.55 , 9.28 , 9.12 , 9.01 , 8.94 , 8.84 , 8.74 , 8.64 , 8.53 } ,
{
7.71 , 6.94 , 6.59 , 6.39 , 6.26 , 6.16 , 6.04 , 5.91 , 5.77 , 5.63 } ,
{
6.61 , 5.79 , 5.41 , 5.19 , 5.05 , 4.95 , 4.82 , 4.68 , 4.53 , 4.36 } ,
{
5.99 , 5.14 , 4.76 , 4.53 , 4.39 , 4.28 , 4.15 , 4 , 3.84 , 3.67 } ,
{
5.59 , 4.74 , 4.35 , 4.12 , 3.97 , 3.87 , 3.73 , 3.57 , 3.41 , 3.23 } ,
{
5.32 , 4.46 , 4.07 , 3.84 , 3.69 , 3.58 , 3.44 , 3.28 , 3.12 , 2.93 } ,
{
5.12 , 4.26 , 3.86 , 3.63 , 3.48 , 3.37 , 3.23 , 3.07 , 2.9 , 2.71 } ,
{
4.96 , 4.1 , 3.71 , 3.48 , 3.33 , 3.22 , 3.07 , 2.91 , 2.74 , 2.54 } ,
{
4.84 , 3.98 , 3.59 , 3.36 , 3.2 , 3.09 , 2.95 , 2.79 , 2.61 , 2.4 } ,
{
4.75 , 3.88 , 3.49 , 3.26 , 3.11 , 3 , 2.85 , 2.69 , 2.5 , 2.3 } ,
{
4.67 , 3.8 , 3.41 , 3.18 , 3.02 , 2.92 , 2.77 , 2.6 , 2.42 , 2.21 } ,
{
4.6 , 3.74 , 3.34 , 3.11 , 2.96 , 2.85 , 2.7 , 2.53 , 2.35 , 2.13 } ,
{
4.54 , 3.68 , 3.29 , 3.06 , 2.9 , 2.79 , 2.64 , 2.48 , 2.29 , 2.07 } ,
{
4.49 , 3.63 , 3.24 , 3.01 , 2.85 , 2.74 , 2.59 , 2.42 , 2.24 , 2.01 } ,
{
4.45 , 3.59 , 3.2 , 2.96 , 2.81 , 2.7 , 2.55 , 2.38 , 2.19 , 1.96 } ,
{
4.41 , 3.55 , 3.16 , 2.93 , 2.77 , 2.66 , 2.51 , 2.34 , 2.15 , 1.92 } ,
{
4.38 , 3.52 , 3.13 , 2.9 , 2.74 , 2.63 , 2.48 , 2.31 , 2.11 , 1.88 } ,
{
4.35 , 3.49 , 3.1 , 2.87 , 2.71 , 2.6 , 2.45 , 2.28 , 2.08 , 1.84 } ,
{
4.32 , 3.47 , 3.07 , 2.84 , 2.68 , 2.57 , 2.42 , 2.25 , 2.05 , 1.81 } ,
{
4.3 , 3.44 , 3.05 , 2.82 , 2.66 , 2.55 , 2.4 , 2.23 , 2.03 , 1.78 } ,
{
4.28 , 3.42 , 3.03 , 2.8 , 2.64 , 2.53 , 2.38 , 2.2 , 2 , 1.76 } ,
{
4.26 , 3.4 , 3.01 , 2.78 , 2.62 , 2.51 , 2.36 , 2.18 , 1.98 , 1.73 } ,
{
4.24 , 3.38 , 2.99 , 2.76 , 2.6 , 2.49 , 2.34 , 2.16 , 1.96 , 1.71 } ,
{
4.22 , 3.37 , 2.98 , 2.74 , 2.59 , 2.47 , 2.32 , 2.15 , 1.95 , 1.69 } ,
{
4.21 , 3.35 , 2.96 , 2.73 , 2.57 , 2.46 , 2.3 , 2.13 , 1.93 , 1.67 } ,
{
4.2 , 3.34 , 2.95 , 2.71 , 2.56 , 2.44 , 2.29 , 2.12 , 1.91 , 1.65 } ,
{
4.18 , 3.33 , 2.93 , 2.7 , 2.54 , 2.43 , 2.28 , 2.1 , 1.9 , 1.64 } ,
{
4.17 , 3.32 , 2.92 , 2.69 , 2.53 , 2.42 , 2.27 , 2.09 , 1.89 , 1.62 } ,
{
4.08 , 3.23 , 2.84 , 2.61 , 2.45 , 2.34 , 2.18 , 2 , 1.79 , 1.51 } ,
{
4 , 3.15 , 2.76 , 2.52 , 2.37 , 2.25 , 2.1 , 1.92 , 1.7 , 1.39 } ,
{
3.92 , 3.07 , 2.68 , 2.45 , 2.29 , 2.17 , 2.02 , 1.83 , 1.61 , 1.25 } ,
{
3.84 , 2.99 , 2.6 , 2.37 , 2.21 , 2.09 , 1.94 , 1.75 , 1.52 , 1 }
};
//α=0.10
static public double[,] F_Dist_10E_2 = {
{
39.86 , 49.5 , 53.59 , 55.83 , 57.24 , 58.2 , 59.44 , 60.71 , 62 , 63.33 } ,
{
8.53 , 9 , 9.16 , 9.24 , 9.29 , 9.33 , 9.37 , 9.41 , 9.45 , 9.49 } ,
{
5.54 , 5.46 , 5.36 , 5.32 , 5.31 , 5.28 , 5.25 , 5.22 , 5.18 , 5.13 } ,
{
4.54 , 4.32 , 4.19 , 4.11 , 4.05 , 4.01 , 3.95 , 3.9 , 3.83 , 3.76 } ,
{
4.06 , 3.78 , 3.62 , 3.52 , 3.45 , 3.4 , 3.34 , 3.27 , 3.19 , 3.1 } ,
{
3.78 , 3.46 , 3.29 , 3.18 , 3.11 , 3.05 , 2.98 , 2.9 , 2.82 , 2.72 } ,
{
3.59 , 3.26 , 3.07 , 2.96 , 2.88 , 2.83 , 2.75 , 2.67 , 2.58 , 2.47 } ,
{
3.46 , 3.11 , 2.92 , 2.81 , 2.73 , 2.67 , 2.59 , 2.5 , 2.4 , 2.29 } ,
{
3.36 , 3.01 , 2.81 , 2.69 , 2.61 , 2.55 , 2.47 , 2.38 , 2.28 , 2.16 } ,
{
3.29 , 2.92 , 2.73 , 2.61 , 2.52 , 2.46 , 2.38 , 2.28 , 2.18 , 2.06 } ,
{
3.23 , 2.86 , 2.66 , 2.54 , 2.45 , 2.39 , 2.3 , 2.21 , 2.1 , 1.97 } ,
{
3.18 , 2.81 , 2.61 , 2.48 , 2.39 , 2.33 , 2.24 , 2.15 , 2.04 , 1.9 } ,
{
3.14 , 2.76 , 2.56 , 2.43 , 2.35 , 2.28 , 2.2 , 2.1 , 1.98 , 1.85 } ,
{
3.1 , 2.73 , 2.52 , 2.39 , 2.31 , 2.24 , 2.15 , 2.05 , 1.94 , 1.8 } ,
{
3.07 , 2.7 , 2.49 , 2.36 , 2.27 , 2.21 , 2.12 , 2.02 , 1.9 , 1.76 } ,
{
3.05 , 2.67 , 2.46 , 2.33 , 2.24 , 2.18 , 2.09 , 1.99 , 1.87 , 1.72 } ,
{
3.03 , 2.64 , 2.44 , 2.31 , 2.22 , 2.15 , 2.06 , 1.96 , 1.84 , 1.69 } ,
{
3.01 , 2.62 , 2.42 , 2.29 , 2.2 , 2.13 , 2.04 , 1.93 , 1.81 , 1.66 } ,
{
2.99 , 2.61 , 2.4 , 2.27 , 2.18 , 2.11 , 2.02 , 1.91 , 1.79 , 1.63 } ,
{
2.97 , 2.59 , 2.38 , 2.25 , 2.16 , 2.09 , 2 , 1.89 , 1.77 , 1.61 } ,
{
2.96 , 2.57 , 2.36 , 2.23 , 2.14 , 2.08 , 1.98 , 1.87 , 1.75 , 1.59 } ,
{
2.95 , 2.56 , 2.35 , 2.22 , 2.13 , 2.06 , 1.97 , 1.86 , 1.73 , 1.57 } ,
{
2.94 , 2.55 , 2.34 , 2.21 , 2.11 , 2.05 , 1.95 , 1.84 , 1.72 , 1.55 } ,
{
2.93 , 2.54 , 2.33 , 2.19 , 2.1 , 2.04 , 1.94 , 1.83 , 1.7 , 1.53 } ,
{
2.92 , 2.53 , 2.32 , 2.18 , 2.09 , 2.02 , 1.93 , 1.82 , 1.69 , 1.52 } ,
{
2.91 , 2.52 , 2.31 , 2.17 , 2.08 , 2.01 , 1.92 , 1.81 , 1.68 , 1.5 } ,
{
2.9 , 2.51 , 2.3 , 2.17 , 2.07 , 2 , 1.91 , 1.8 , 1.67 , 1.49 } ,
{
2.89 , 2.5 , 2.29 , 2.16 , 2.06 , 2 , 1.9 , 1.79 , 1.66 , 1.48 } ,
{
2.89 , 2.5 , 2.28 , 2.15 , 2.06 , 1.99 , 1.89 , 1.78 , 1.65 , 1.47 } ,
{
2.88 , 2.49 , 2.28 , 2.14 , 2.05 , 1.98 , 1.88 , 1.77 , 1.64 , 1.46 } ,
{
2.84 , 2.44 , 2.23 , 2.09 , 2 , 1.93 , 1.83 , 1.71 , 1.57 , 1.38 } ,
{
2.79 , 2.39 , 2.18 , 2.04 , 1.95 , 1.87 , 1.77 , 1.66 , 1.51 , 1.29 } ,
{
2.75 , 2.35 , 2.13 , 1.99 , 1.9 , 1.82 , 1.72 , 1.6 , 1.45 , 1.19 } ,
{
2.71 , 2.3 , 2.08 , 1.94 , 1.85 , 1.17 , 1.67 , 1.55 , 1.38 , 1 } ,
};
/// <summary>
/// 常规最小二乘计算回归系数、常数
/// </summary>
/// <param name="x"></param>
/// <param name="y"></param>
/// <returns></returns>
public static double[,] OLSMulLinearRegression(double[,] x, double[,] y)
{
double[,] xTrans = Matrix.Transpose(x);
double[,] xTx = Matrix.MultiplyMatrix(xTrans, x);
double[,] xTxAth = Matrix.Athwart(xTx);
double[,] beta = Matrix.MultiplyMatrix(xTxAth, xTrans);
beta = Matrix.MultiplyMatrix(beta, y);
return beta;
}
/// <summary>
/// 偏F显著性统计量
/// 此处x样例首列为1
/// </summary>
/// <param name="x">此处x样例首列为1,已备归常量</param>
/// <param name="y">因变量元数据</param>
/// <returns></returns>
public static double[] PartialFSigStatistic(double[,] x, double[,] y)
{
int n = y.GetLength(0);
int p = x.GetLength(1) - 1;//此处为添加首列1后x矩阵
double[] retval = new double[p];//返回偏F统计量
double[] t_Statistics = new double[p];
double[,] beta = OLSMulLinearRegression(x, y);
double[,] yPredict = Matrix.MultiplyMatrix(x, beta);
double yAvg = Matrix.SubDatas(y, 0, RowGet: false).Average();
double sse = 0;
double sst = 0;
for (int i = 0; i < n; i++)
{
sst += Math.Pow(y[i, 0] - yAvg, 2);
sse += Math.Pow(y[i, 0] - yPredict[i, 0], 2);
}
for (int i = 1; i < p + 1; i++)
{
double[,] ix = Matrix.MatrixRemoveVector(x, i);
double[,] iBeta = OLSMulLinearRegression(ix, y);
double[,] iPreY = Matrix.MultiplyMatrix(ix, iBeta);
double isse = 0;
for (int j = 0; j < n; j++)
{
isse += Math.Pow(y[j, 0] - iPreY[j, 0], 2);
}
retval[i - 1] = (isse - sse) / (sse / (n - p - 1));
t_Statistics[i - 1] = Math.Sqrt(retval[i - 1]);
}
double F_Statistics = ((sst - sse) / p) / (sse / (n - p - 1));
return retval;
}
/// <summary>
/// 样本数据标准化
/// x[r,c]=(x[r,c]-avg(x[c]))/(Σ(x[r,c]-avg(x[c]))^(1/2))
/// </summary>
/// <param name="MatrixEin">样本矩阵</param>
/// <returns></returns>
public static double[,] SampleDataStandardized(double[,] MatrixEin)
{
int eRow = MatrixEin.GetLength(0);
int eCol = MatrixEin.GetLength(1);
double[,] outMatrix = new double[eRow, eCol];
double[] colAverages = new double[eCol];
double[] colSqrErrsSum = new double[eCol];
for (int i = 0; i < eCol; i++)
{
colAverages[i] = Matrix.SubDatas(MatrixEin, i, RowGet: false).Average();
}
for (int i = 0; i < eRow; i++)
{
for (int j = 0; j < eCol; j++)
{
colSqrErrsSum[j] += Math.Pow(MatrixEin[i, j] - colAverages[j], 2);
}
}
for (int r = 0; r < eRow; r++)
{
for (int c = 0; c < eCol; c++)
{
outMatrix[r, c] = (MatrixEin[r, c] - colAverages[c]) / Math.Sqrt(colSqrErrsSum[c]);
}
}
return outMatrix;
}
/// <summary>
/// 样本相关系数矩阵
/// r=(X*T)(X*)
/// </summary>
/// <param name="MatrixEin">标准化后样本数据</param>
/// <returns></returns>
public static double[,] SampleRelativeCoeffMatrix(double[,] MatrixEin)
{
double[,] traMatrixEin = Matrix.Transpose(MatrixEin);
double[,] outMatrix = Matrix.MultiplyMatrix(traMatrixEin, MatrixEin);
return outMatrix;
}
/// <summary>
/// 逐步回归
/// 可用于有效自变量选择,
/// 多重共线性未解决
/// </summary>
/// <param name="x">已添加首列为1自变量矩阵</param>
/// <param name="y"></param>
/// <returns>有效自变量列id</returns>
public static int[] StepwiseRegression(double[,] x, double[,] y)
{
int n = y.GetLength(0);
int p = x.GetLength(1) - 1;
//获取除去首列(全为1)的列,方便标准化样本
double[,] noneX = Matrix.MatrixRemoveVector(x, 0, RowRmv: false);
//x,y样本组合阵
double[,] xyGroup = Matrix.MatrixCombine(y, noneX);
//xy组合阵标准化
double[,] xyStandGroup = SampleDataStandardized(xyGroup);
//xy样本相关系数阵
double[,] xyRelaCoeffs = SampleRelativeCoeffMatrix(xyStandGroup);
//xy样本相关系数首行
double[] xy1Row = Matrix.SubDatas(xyRelaCoeffs, 0, RowGet: true);
//xy样本首行降序排序(当前首元素为1,不需F统计量偏检测)
double[] xy1RowDecs = xy1Row.ToList().OrderByDescending(x => x).ToArray();
//依据自变量与因变量之间的相关系数
//从大到小分别插入回归公式
//进行偏F统计量检测:
// Fentry>=Fs.t.entray-引入,否则不引入,此次引入α-5%
// Fremoval<=Fs.t.removal-剔除,否则不剔除,此次剔除α-10%
double[,] validData = Matrix.SubMatrix(x, 0, RowGet: false);
List<int> validColId = new List<int>();
//int[] testSeq = { 0, 2, 6, 1, 4, 3, 5 };
for (int i = 1; i < xy1RowDecs.Length; i++)
{
//获取通过相关系数排序后列索引(×此处+1:该处使用xy1Row为无首列1矩阵)
int colId = xy1Row.ToList().IndexOf(xy1RowDecs[i]);
//int colId = testSeq[i];
validColId.Add(colId);
//获取该id处的元数据
double[] colDd = Matrix.SubDatas(x, colId, RowGet: false);
//将选择列数据与有效数据合并
validData = Matrix.MatrixCombine(validData, Matrix.VectorGenerate(colDd, RowVec: false));
//计算当前有效数据与y元数据偏F统计量
double[] partFStatistics = PartialFSigStatistic(validData, y);
//获取当前有效数据中自变量数量
int currP = validData.GetLength(1) - 1;
//F表查询-引入临界F值α-1%
double f_Entry = F_Dist_5E_2[n - currP - 2, currP - 1];
//F表查询-剔除临界F值α-2.5%
double f_Removal = F_Dist_10E_2[n - currP - 2, currP - 1];
//如果该自变量计算所得偏F统计量不满足引入临界F值则不引入
if (partFStatistics[partFStatistics.Length - 1] < f_Entry)
{
validData = Matrix.MatrixRemoveVector(validData, validData.GetLength(1) - 1, RowRmv: false);
validColId.Remove(colId);
continue;
}
//如果该自变量计算所得偏F统计量满足引入临界F值时,对当前全部自变量进行剔除筛选
//待剔除列id,在此基于临时有效数据
List<int> removalIds = new List<int>();
//从新引入后的有效数据中筛选待剔除和无需引入的数据
for (int j = 0; j < currP; j++)
{
//当j项自变量偏F统计量满足剔除条件时
if (partFStatistics[j] <= f_Removal)
{
//记录待剔除项列在有效数据中的索引
removalIds.Add(j);
//删除之前记录的有效数据对应元数据列索引
validColId.RemoveAt(j);
}
}
if (removalIds.Count != 0)
{
removalIds = removalIds.OrderByDescending(x => x).ToList();
for (int j = 0; j < removalIds.Count; j++)
{
//注:当前待剔除索引基准时自变量基准,
//因有效数据中自变量基于1开始,因此此处剔除时需要索引+1
validData = Matrix.MatrixRemoveVector(validData, removalIds[j] + 1, RowRmv: false);
}
}
}
return validColId.ToArray();
}
/// <summary>
/// 决定系数R²=SSR/SST=1-SSE/SST
/// </summary>
/// <param name="x"></param>
/// <param name="y"></param>
/// <returns></returns>
public static double DecisionCoefficient(double[,] x, double[,] y)
{
int n = y.GetLength(0);
double yAvg = Matrix.SubDatas(y, 0, RowGet: false).ToList().Average();
double[,] beta = OLSMulLinearRegression(x, y);
double[,] yPre = Matrix.MultiplyMatrix(x, beta);
double SST = 0;
double SSR = 0;
for (int i = 0; i < n; i++)
{
SST += Math.Pow(y[i, 0] - yAvg, 2);
SSR += Math.Pow(yPre[i, 0] - yAvg, 2);
}
return (double)SSR / SST;
}
/// <summary>
/// 样本偏决定系数
/// 公式:R²(j)=[SSE(j)-SSE]/SSE,
/// 偏相关系数:R(j)=(R²)^(1/2)
/// </summary>
/// <param name="x"></param>
/// <param name="y"></param>
/// <returns></returns>
public static double[] PartialDecisionCoefficient(double[,] x, double[,] y)
{
int p = x.GetLength(1);
int n = y.GetLength(0);
double[] partDecCoes = new double[p];
double[,] beta = OLSMulLinearRegression(x, y);
double[,] yPredict = Matrix.MultiplyMatrix(x, beta);
double SSE = 0;
for (int i = 0; i < n; i++)
{
SSE += Math.Pow(y[i, 0] - yPredict[i, 0], 2);
}
for (int i = 0; i < p; i++)
{
double[,] ix = Matrix.MatrixRemoveVector(x, i, RowRmv: false);
double[,] iyPre = Matrix.MultiplyMatrix(ix, OLSMulLinearRegression(ix, y));
double isse = 0;
for (int j = 0; j < n; j++)
{
isse += Math.Pow(y[j, 0] - iyPre[j, 0], 2);
}
partDecCoes[i] = (isse - SSE) / isse;
}
return partDecCoes;
}
/// <summary>
/// 方差扩大因子:
/// 可用于判定解释变量样本多重共线性,生成矩阵对角线元素为对应列的方差扩大因子;
/// 公式:VIF=(X*TX*)^-1;
/// 要求:输入解释变量矩阵首列不为1
/// </summary>
/// <param name="x"></param>
/// <returns></returns>
public static double[,] VarianceInflationFactor(double[,] x)
{
double[,] standX = SampleDataStandardized(x);
double[,] standXT = Matrix.Transpose(standX);
double[,] VIF = Matrix.MultiplyMatrix(standXT, standX);
return Matrix.Athwart(VIF);
}
/// <summary>
/// ret=1/(1+e^(-x))
/// </summary>
/// <param name="x"></param>
/// <returns></returns>
public static double Sigmoid(double x)
{
double[,] a = {
{
1 }, {
1 } };
return 1.0 / (1 + Math.Exp(-1.0 * x));
}
/// <summary>
/// 计算逻辑斯特梯度
/// </summary>
/// <param name="x"></param>
/// <param name="y"></param>
/// <param name="beta"></param>
/// <returns></returns>
public static double[,] LogisticGradientCalculate(double[,] x, double[,] y, double[,] beta)
{
int n = x.GetLength(0);
int p = x.GetLength(1);
double[,] err = new double[n, 1];
for (int i = 0; i < n; i++)
{
double[,] currRow = Matrix.SubMatrix(x, i, RowGet: true);
err[i, 0] = y[i, 0] - Sigmoid(Matrix.MultiplyMatrix(currRow, beta)[0, 0]);
}
double[,] xTra = Matrix.Transpose(x);
double[,] dP_dB = Matrix.MultiplyMatrix(xTra, err);
return Matrix.MultiplyConst(dP_dB, 1.0 / (double)n);
}
/// <summary>
/// 梯度下降计算逻辑回归系数
/// </summary>
/// <param name="x"></param>
/// <param name="y"></param>
/// <param name="alpha">学习率</param>
/// <param name="iterations">迭代次数</param>
/// <returns></returns>
public static double[,] GradientDescentInLogisticRegression(double[,] x, double[,] y, double alpha, int iterations)
{
int p = x.GetLength(1);
double[,] beta = new double[p, 1];
for (int i = 0; i < iterations; i++)
{
//当前β计算logistic相应损失函数梯度(各参数偏导)
double[,] dl = LogisticGradientCalculate(x, y, beta);
//通过每次递减,找寻β系数迭代递减结果
beta = Matrix.AddMatrix(beta, Matrix.MultiplyConst(dl, alpha));
}
return beta;
}
public static double[,] FGDInLogisticRegression(double[,] x, double[,] y, double alpha, int iterations, double epsilon = 1e-6)
{
int n = x.GetLength(0);
int p = x.GetLength(1);
double[,] xTran = Matrix.Transpose(x);
double[,] beta = new double[p, 1];
double norm = 0;
for (int t = 0; t < iterations; t++)
{
double[,] error = new double[n, 1];
//计算当前beta系数下,损失函数对未知参数梯度(在此选择全部样本)
for (int i = 0; i < n; i++)
{
double[,] currRow = Matrix.SubMatrix(x, i, RowGet: true);
error[i, 0] = y[i, 0] - Sigmoid(Matrix.MultiplyMatrix(currRow, beta)[0, 0]);
}
double[,] dl = new double[p, 1];
dl = Matrix.MultiplyConst(Matrix.MultiplyMatrix(xTran, error), 1.0 / (double)n);//梯度向量
norm = Matrix.Norm(dl, NormType.L1);
if (norm < epsilon)
{
break;
}
beta = Matrix.AddMatrix(beta, Matrix.MultiplyConst(dl, alpha));//梯度递减迭代
}
return beta;
}
/// <summary>
/// 随机梯度下降感知机
/// </summary>
/// <param name="x">样本矩阵(存在首列1向量)</param>
/// <param name="y">标签矩阵(列)</param>
/// <param name="lr">学习率</param>
/// <param name="iterations">最大迭代次数</param>
/// <returns></returns>
public static double[,] PerceptronModel(double[,] x, double[,] y, double lr = 1e-2, double iterations = 1e4)
{
int p = x.GetLength(1);
int n = x.GetLength(0);
double[,] beta = new double[p, 1];
//初始化权重矩阵[0,0,…,0]
int rId = 0;
for (int t = 0; t < iterations; t++)
{
rId = new Random().Next(n);
//随机找寻样本数据相应索引
double[,] currRow = Matrix.SubMatrix(x, rId, RowGet: true);
//获取随机索引样本数据
double linePart = Matrix.MultiplyMatrix(currRow, beta)[0, 0];
//理论超平面值计算
bool check = -1 * y[rId, 0] * Sign(linePart) >= 0;
//误分类判断
if (!check)
{
continue;
}
else
{
double[,] lossGrad = Matrix.MultiplyConst(currRow, y[rId, 0]);
//损失函数梯度计算
lossGrad = Matrix.MultiplyConst(lossGrad, lr);
//学习率*损失函数梯度
beta = Matrix.AddMatrix(beta, Matrix.Transpose(lossGrad));
//权重系数更新
}
}
return beta;
}
/// <summary>
/// Sign函数
/// </summary>
/// <param name="x"></param>
/// <returns>x > 0 ? 1 : x == 0 ? 0 : -1</returns>
public static int Sign(double x)
{
return x > 0 ? 1 : x == 0 ? 0 : -1;
}
Func<double[,], double[,], double[,], double, double, double[,]> AdaGrad = (weight, S, grad, eta, epsilon) =>
{
return S;
};
}
}
c#回归相关算法
猜你喜欢
转载自blog.csdn.net/JAYLEE900/article/details/132508133
今日推荐
周排行