C# 拟合曲线项目：PolyFitCSharp
C# 拟合曲线参考资料：
using System; using System.Collections.Generic; using System.Collections.ObjectModel; using System.Linq; using System.Text; using System.Threading.Tasks; using System.Windows.Controls; using System.Windows.Data; using System.Windows.Documents; using System.Windows.Input; using System.Windows.Media; using System.Windows.Media.Imaging; using System.Windows.Navigation; using System.Windows.Shapes; using Microsoft.Research.DynamicDataDisplay; using Microsoft.Research.DynamicDataDisplay.DataSources; using Microsoft.Research.DynamicDataDisplay.PointMarkers; using ViliPetek.LinearAlgebra; using MathNet; using MathNet.Numerics; using MathNet.Numerics.Differentiation; using MathNet.Numerics.Distributions; using MathNet.Numerics.Financial; using MathNet.Numerics.IntegralTransforms; using MathNet.Numerics.Integration; using MathNet.Numerics.Interpolation; using MathNet.Numerics.LinearAlgebra; using MathNet.Numerics.LinearRegression; using MathNet.Numerics.OdeSolvers; using MathNet.Numerics.Optimization; using MathNet.Numerics.Properties; using MathNet.Numerics.Providers; using MathNet.Numerics.Random; using MathNet.Numerics.RootFinding; using MathNet.Numerics.Statistics; using MathNet.Numerics.LinearAlgebra.Double; using MathNet.Numerics.LinearAlgebra.Factorization; namespace PolyFitExample { ///// <summary> ///// 摘自：https://stackoverflow.com/questions/7792088/non-linear-regression-in-c-sharp ///// </summary> //public class PolynomialRegression1 //{ // Vector x_data, y_data, coef; // int order; // public PolynomialRegression1(Vector x_data, Vector y_data, int order) // { // if (x_data.Count != y_data.Count) // { // throw new IndexOutOfRangeException(); // } // this.x_data = x_data; // this.y_data = y_data; // this.order = order; // int N = x_data.Count; // MathNet.Numerics.LinearAlgebra.Double.Matrix A = new MathNet.Numerics.LinearAlgebra.Double.Matrix(N, order + 1); // for (int i = 0; i < N; i++) // { // A.SetRow(VandermondeRow(x_data[i]), i); // } // // Least Squares of |y=A(x)*c| // // tr(A)*y = tr(A)*A*c // // inv(tr(A)*A)*tr(A)*y = c // Matrix At = Matrix.Transpose(A); // Matrix y2 = new Matrix(y_data, N); // coef = (At * A).Solve(At * y2).GetColumnVector(0); // } // Vector VandermondeRow(double x) // { // double[] row = new double[order + 1]; // for (int i = 0; i <= order; i++) // { // row[i] = Math.Pow(x, i); // } // return new DenseVector(row); // } // public double Fit(double x) // { // return Vector.ScalarProduct(VandermondeRow(x), coef); // } // public int Order { get { return order; } } // public Vector Coefficients { get { return coef; } } // public Vector XData { get { return x_data; } } // public Vector YData { get { return y_data; } } //} }
C# 拟合曲线实际应用：
using System; using System.Collections.Generic; using System.Collections.ObjectModel; using System.Linq; using System.Text; using System.Threading.Tasks; using System.Windows; using System.Windows.Controls; using System.Windows.Data; using System.Windows.Documents; using System.Windows.Input; using System.Windows.Media; using System.Windows.Media.Imaging; using System.Windows.Navigation; using System.Windows.Shapes; using Microsoft.Research.DynamicDataDisplay; using Microsoft.Research.DynamicDataDisplay.DataSources; using Microsoft.Research.DynamicDataDisplay.PointMarkers; using ViliPetek.LinearAlgebra; using MathNet; using MathNet.Numerics; using MathNet.Numerics.Differentiation; using MathNet.Numerics.Distributions; using MathNet.Numerics.Financial; using MathNet.Numerics.IntegralTransforms; using MathNet.Numerics.Integration; using MathNet.Numerics.Interpolation; using MathNet.Numerics.LinearAlgebra; using MathNet.Numerics.LinearRegression; using MathNet.Numerics.OdeSolvers; using MathNet.Numerics.Optimization; using MathNet.Numerics.Properties; using MathNet.Numerics.Providers; using MathNet.Numerics.Random; using MathNet.Numerics.RootFinding; using MathNet.Numerics.Statistics; using MathNet.Numerics.LinearAlgebra.Double; using MathNet.Numerics.LinearAlgebra.Factorization; namespace PolyFitExample { /// <summary> /// Interaction logic for MainWindow.xaml /// </summary> public partial class MainWindow : System.Windows.Window { private ObservableCollection<Point> _dataSeries; private ObservableCollection<Point> _fittedSeries; public MainWindow() { InitializeComponent(); InitializeGraphs(); } private void InitializeGraphs() { // create the storage for the series _dataSeries = new ObservableCollection<Point>(); _fittedSeries = new ObservableCollection<Point>(); // create the series chartPlotter.AddLineGraph(_dataSeries.AsDataSource(), new Pen(Brushes.Blue, 1), new CircleElementPointMarker { Size = 3, Brush = Brushes.Blue, Fill = Brushes.Blue }, new PenDescription("Data")); chartPlotter.AddLineGraph(_fittedSeries.AsDataSource(), new Pen(Brushes.Red, 3), new PenDescription("PolyFit")); // hide the standard legend chartPlotter.LegendVisibility = System.Windows.Visibility.Hidden; chartPlotter.NewLegendVisible = false; } private void graphButton_Click(object sender, RoutedEventArgs e) { const int DataSize = 1024; double[] x = new double[DataSize]; double[] y = new double[DataSize]; Random rand = new Random(); for (int i = 0; i < DataSize; i++) { x[i] = i; y[i] = Math.Sin(Math.PI * (i / (double)180)); //Math.Sin((double)i / 1024.0 * Math.PI * 2) + (rand.NextDouble() - 0.5); } //方法一 //Vector<double> pp = Polyfit01(x, y, (int)orderField.Value); //List<double> fity = new List<double>(); //foreach (int xitem in x) //{ // double yitem = Polyval(pp, xitem, (int)orderField.Value); // fity.Add(yitem); //} //方法二 //Vector<double> desRes = Polyfit02(x, y, (int)orderField.Value); //List<double> fity = new List<double>(); //foreach (int xitem in x) //{ // double yitem = Polyval(desRes, xitem, (int)orderField.Value); // fity.Add(yitem); //} //方法三 Vector<double> desRes = Polyfit03(x, y, (int)orderField.Value); List<double> fity = new List<double>(); foreach (int xitem in x) { double yitem = Polyval(desRes, xitem, (int)orderField.Value); fity.Add(yitem); } // plot the data _dataSeries.SetXY(x, y); _fittedSeries.SetXY(x, fity.ToArray()); } public static Vector<double> Polyfit01(DenseVector xData, DenseVector yData, int order) { DenseMatrix vandMatrix = new DenseMatrix(xData.Count, order + 1); for (int i = 0; i < xData.Count; i++) { double[] rowResult = new double[vandMatrix.ColumnCount]; double mult = 1; for (int rowIndex = 0; rowIndex < rowResult.Count(); rowIndex++) { rowResult[rowIndex] = mult; mult *= xData[i]; } vandMatrix.SetRow(i, new DenseVector(rowResult)); } DenseVector result = new DenseVector(vandMatrix.ColumnCount); for (int j = 0; j < vandMatrix.RowCount; j++) { for (int i = 0; i < vandMatrix.ColumnCount; i++) { result[i] += vandMatrix[j, i] * yData[j]; } } return vandMatrix.TransposeThisAndMultiply(vandMatrix).LU().Solve(result); } public static Vector<double> Polyfit02(double[] x, double[] y, int order) { Matrix<double> design = Matrix<double>.Build.Dense(x.Length, order + 1, (i, j) => Math.Pow(x[i], j)); return MultipleRegression.QR(design, Vector<double>.Build.Dense(y)); } public static Vector<double> Polyfit03(double[] x, double[] y, int degree) { DenseMatrix v = new DenseMatrix(x.Length, degree + 1); for (int i = 0; i < v.RowCount; i++) for (int j = 0; j <= degree; j++) v[i, j] = Math.Pow(x[i], j); Matrix<double> yv = new DenseVector(y).ToColumnMatrix(); QR<double> qr = v.QR(); var r = qr.R.SubMatrix(0, degree + 1, 0, degree + 1); var q = v.Multiply(r.Inverse()); var p = r.Inverse().Multiply(q.TransposeThisAndMultiply(yv)); return p.Column(0); } public static double Polyval(Vector<double> _coefs, double x, int order) { double[] result = new double[order + 1]; double mult = 1; for (int i = 0; i < result.Count(); i++) { result[i] = mult; mult *= x; } return _coefs * new DenseVector(result); } //polyfit用于多项式曲线拟合 //p = polyfit(x, y, m) //其中, x, y为已知数据点向量, 分别表示横,纵坐标, m为拟合多项式的次数, 结果返回m次拟合多项式系数, 从高次到低次存放在向量p中. //y0=polyval(p, x0) //可求得多项式在x0处的值y0 //public static void TestGet() //{ // const int degree = 2; // var x = new[] { 0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0 }; // var y = new[] { 1.0, 6.0, 17.0, 34.0, 57.0, 86.0, 121.0, 162.0, 209.0, 262.0, 321.0 }; // var p = Polyfit(x, y, degree); // foreach (var d in p) Console.Write("{0} ", d); // Console.WriteLine(); // for (int i = 0; i < x.Length; i++) // Console.WriteLine("{0} => {1} diff {2}", x[i], Polyval(p, x[i]), y[i] - Polyval(p, x[i])); // Console.ReadKey(true); //} } }