使用MathNet.Numerics进行曲线拟合并使用SkiaSharp绘制曲线

本文展示了如何使用MathNet.Numerics对离散点进行曲线拟合，计算其R^2值，并基于Winform使用SkiaSharp.Views.WindowsForms绘制曲线及离散点，上述组件可从NuGet获取。

1、曲线拟合

使用Fit.Polynomial方法获取拟合后的多项式系数数组，其参数分别为x坐标数组、y坐标数组、曲线阶数，返回值用于计算多项式的值。

//曲线拟合
var x = _skOriginalPointList.Select(t => (double)t.X).ToArray();
var y = _skOriginalPointList.Select(t => (double)t.Y).ToArray();
_order = Math.Min(5, x.Length - 1); //曲线阶数
var ret = Fit.Polynomial(x, y, _order);

 2、多项式的值计算

多项式的值计算方法：p0 + p1*x + p2*x^2 + ... + pk*x^k。

/// <summary>
/// 获取多项式的值
/// </summary>
/// <param name="order">阶数</param>
/// <param name="x">x坐标</param>
/// <param name="polynomialArray">多项式系数数组</param>
/// <returns>Y坐标</returns>
private static double GetPolynomialValue(int order, double x, double[] polynomialArray) =>
    Enumerable.Range(0, order + 1)
        .Select(t => polynomialArray[t] * Math.Pow(x, t))
        .Sum();

3、获取曲线

曲线绘制的宽高分别为width、height，那么x坐标从0开始，每间隔一定像素获取一个点，作为绘制曲线的点集合。

//获取拟合后的曲线
_skFittingPointList.Clear();
for (var i = 0; i <= width; i++)
{
    var value = GetPolynomialValue(_order, i, ret);
    _skFittingPointList.Add(new SKPoint(i, (float)value));
    i += 5;
}

4、计算R^2值

//计算R^2,R^2这个值越接近1，说明拟合出来的曲线跟原曲线就越接近
var yTest = _skOriginalPointList.Select(point =>
        GetPolynomialValue(_order, point.X, ret))
    .ToArray();
_rSquared = GoodnessOfFit.RSquared(y, yTest);

5、曲线图展示

6、主要代码 

using System;
using System.Collections.Generic;
using System.Drawing;
using System.Linq;
using System.Windows.Forms;
using MathNet.Numerics;
using SkiaSharp;
using SkiaSharp.Views.Desktop;

namespace MathSample
{
    public partial class FormMain : Form
    {
        private readonly Point[] _originPoints;
        private readonly List<SKPoint> _skOriginalPointList;
        private readonly List<SKPoint> _skFittingPointList;

        private readonly Panel _panel;
        private readonly SKControl _skiaView;

        private int _order; //曲线阶数
        private double _rSquared; //R^2值

        public FormMain()
        {
            InitializeComponent();

            //原始点(百分比)
            _originPoints = new[]
            {
                new Point(16, 60),
                new Point(40, 20),
                new Point(52, 47),
                new Point(63, 85),
                new Point(87, 26)
            };
            _skOriginalPointList = new List<SKPoint>();
            _skFittingPointList = new List<SKPoint>();

            _panel = new Panel();

            //skiaView
            _skiaView = new SKControl
            {
                Dock = DockStyle.Fill,
                Location = new Point(0, 0),
                Name = "skiaView",
                TabIndex = 0
            };
            _skiaView.PaintSurface += skiaView_PaintSurface;
        }

        private void FormMain_Load(object sender, EventArgs e)
        {
            Controls.Add(_panel);
            _panel.Controls.Add(_skiaView);
            ResizePanel();
            DrawLines();
        }

        private void FormMain_SizeChanged(object sender, EventArgs e)
        {
            ResizePanel();
            DrawLines();
        }

        private void ResizePanel()
        {
            _panel.Left = 10;
            _panel.Top = 10;
            _panel.Width = ClientSize.Width - 20;
            _panel.Height = ClientSize.Height - 20;
        }

        private void DrawLines()
        {
            var width = _panel.Width;
            var height = _panel.Height;
            _skOriginalPointList.Clear();
            _skOriginalPointList.AddRange(_originPoints.Select(t =>
                new SKPoint(t.X / 100f * width, (100 - t.Y) / 100f * height)));

            //曲线拟合
            var x = _skOriginalPointList.Select(t => (double)t.X).ToArray();
            var y = _skOriginalPointList.Select(t => (double)t.Y).ToArray();
            _order = Math.Min(5, x.Length - 1); //曲线阶数
            var ret = Fit.Polynomial(x, y, _order);

            //获取拟合后的曲线
            _skFittingPointList.Clear();
            for (var i = 0; i <= width; i++)
            {
                var value = GetPolynomialValue(_order, i, ret);
                _skFittingPointList.Add(new SKPoint(i, (float)value));
                i += 5;
            }

            //计算R^2,R^2这个值越接近1，说明拟合出来的曲线跟原曲线就越接近
            var yTest = _skOriginalPointList.Select(point =>
                    GetPolynomialValue(_order, point.X, ret))
                .ToArray();
            _rSquared = GoodnessOfFit.RSquared(y, yTest);

            //刷新画布
            _skiaView.Refresh();
        }

        /// <summary>
        /// 获取多项式的值
        /// </summary>
        /// <param name="order">阶数</param>
        /// <param name="x">x坐标</param>
        /// <param name="polynomialArray">多项式系数数组</param>
        /// <returns>Y坐标</returns>
        private static double GetPolynomialValue(int order, double x, double[] polynomialArray) =>
            Enumerable.Range(0, order + 1)
                .Select(t => polynomialArray[t] * Math.Pow(x, t))
                .Sum();

        private void skiaView_PaintSurface(object sender, SKPaintSurfaceEventArgs e)
        {
            // the the canvas and properties
            var canvas = e.Surface.Canvas;

            //绘图
            DrawSkia(canvas);
        }

        private void DrawSkia(SKCanvas canvas)
        {
            // get the screen density for scaling
            var scale = 1f;

            // handle the device screen density
            canvas.Scale(scale);

            // make sure the canvas is blank
            canvas.Clear(GetSkColor(Color.LightGray));

            DrawLines(canvas, _skFittingPointList, Color.Aqua);
            DrawPoints(canvas, _skOriginalPointList, Color.Red);
            DrawText(canvas, Color.Black);
        }

        //画曲线
        private void DrawLines(SKCanvas canvas, IList<SKPoint> points, Color color)
        {
            if (points.Count == 0)
            {
                return;
            }

            using (var paint = new SKPaint())
            {
                paint.Style = SKPaintStyle.Stroke;
                paint.StrokeWidth = 3;
                paint.IsAntialias = true; //抗锯齿
                paint.StrokeCap = SKStrokeCap.Round;
                paint.Color = GetSkColor(color);

                using (var path = new SKPath())
                {
                    path.MoveTo(points[0]);
                    for (var i = 1; i < points.Count; i++)
                    {
                        path.LineTo(points[i]);
                    }

                    canvas.DrawPath(path, paint);
                }
            }
        }

        //画点
        private void DrawPoints(SKCanvas canvas, List<SKPoint> pointList, Color color)
        {
            using (var paint = new SKPaint())
            {
                paint.Style = SKPaintStyle.Stroke;
                paint.StrokeWidth = 8;
                paint.IsAntialias = true; //抗锯齿
                paint.StrokeCap = SKStrokeCap.Round;
                paint.Color = GetSkColor(color);
                pointList.ForEach(t => canvas.DrawPoint(t, paint));
            }
        }

        //文本
        private void DrawText(SKCanvas canvas, Color color)
        {
            // draw some text
            var paint = new SKPaint
            {
                Color = GetSkColor(color),
                IsAntialias = true,
                Style = SKPaintStyle.Fill,
                TextAlign = SKTextAlign.Left,
                TextSize = 16
            };
            var point = new SKPoint(5, paint.TextSize + 5);
            var text = $"{_order}阶曲线拟合：R^2 = {_rSquared}";
            canvas.DrawText(text, point, paint);
        }

        private static SKColor GetSkColor(Color color)
        {
            return new SKColor(color.R, color.G, color.B);
        }
    }
}