using MathNet.Numerics;
|
using MathNet.Numerics.LinearAlgebra.Double;
|
using NPOI.SS.Formula.Functions;
|
using System;
|
using System.Collections.Generic;
|
using System.Drawing;
|
using System.Linq;
|
using System.Text;
|
using System.Threading.Tasks;
|
using System.Windows.Forms.DataVisualization.Charting;
|
|
namespace Hydro.CommonBase
|
{
|
public class CurveFitHelper
|
{
|
|
|
public DenseVector coefficients;
|
int degree;
|
double xMax;
|
double xMin;
|
List<PointF> data;
|
bool isFitted = false;
|
|
public CurveFitHelper(List<PointF> data, int degree = 2)
|
{
|
this.degree = degree;
|
this.data = data;
|
|
double[] xlist = data.Select(p => (double)p.X).ToArray();
|
if (data.Count>2)
|
{
|
coefficients = Fit.Polynomial(xlist, data.Select(p => (double)p.Y).ToArray(), 2);
|
isFitted = true;
|
}
|
else if(data.Count==2)
|
{
|
//使用两个点的线性拟合
|
var tuple = Fit.Line(xlist, data.Select(p => (double)p.Y).ToArray());
|
coefficients= new DenseVector(new double[] { tuple.Item1, tuple.Item2 });
|
isFitted = true;
|
}
|
else if (data.Count==1)
|
{
|
coefficients = new DenseVector(new double[] { data[0].Y });
|
isFitted = true;
|
}
|
else
|
{
|
coefficients = new DenseVector(new double[] { 0});
|
isFitted = true;
|
}
|
|
var list = xlist.ToList();
|
xMax = list.Max();
|
xMin = list.Min();
|
}
|
|
public double Evaluate(double x)
|
{
|
if (!isFitted) return double.NaN;
|
double y = 0;
|
for (int i = 0; i < coefficients.Count; i++)
|
{
|
y += (float)(coefficients[i] * Math.Pow(x, i));
|
}
|
return y;
|
}
|
|
|
public List<PointF> GetFitCurve(int count)
|
{
|
List<PointF> curve = new List<PointF>();
|
|
for (double x = xMin; x <= xMax && curve.Count<count; x += (xMax - xMin) / count)
|
{
|
double y = Evaluate(x);
|
curve.Add(new PointF((float)x, (float)y));
|
}
|
return curve;
|
|
}
|
|
public List<List<PointF>> GetConfidenceCurve(int count)
|
{
|
List<PointF> up_list = new List<PointF>();
|
List<PointF> down_list = new List<PointF>();
|
List<PointF> avg_list = new List<PointF>();
|
|
double X =0;
|
double Y=0;
|
for (double x = xMin; x < xMax; x += (xMax - xMin) / count)
|
{
|
var list = data.FindAll(p => x - (xMax - xMin) / count / 2 <= p.X && p.X < x + (xMax - xMin) / count/2).Select(p => (double)p.Y).ToList();
|
if (list.Count >= 3)
|
{
|
var Bounds = GetConfidenceRange(list, 95);
|
var x_avg = x;
|
if (Bounds[0] == double.NaN || Bounds[1] == double.NaN) continue;
|
down_list.Add(new PointF((float)x_avg, (float)Bounds[0]));
|
up_list.Add(new PointF((float)x_avg, (float)Bounds[1]));
|
X = x;
|
//平均值
|
Y = (Bounds[0] + Bounds[1] )/ 2;
|
avg_list.Add(new PointF((float)x_avg, (float)Y));
|
}
|
else if (list.Count == 2)
|
{
|
var x_avg = x;
|
down_list.Add(new PointF((float)x_avg, (float)list.Min()));
|
up_list.Add(new PointF((float)x_avg, (float)list.Max()));
|
X = x;
|
Y = (list.Min() + list.Max()) / 2;
|
avg_list.Add(new PointF((float)x_avg, (float)Y));
|
}
|
else if (list.Count == 1)
|
{
|
var x_avg = x;
|
down_list.Add(new PointF((float)x_avg, (float)list[0]));
|
up_list.Add(new PointF((float)x_avg, (float)list[0]));
|
X = x;
|
Y = list[0];
|
avg_list.Add(new PointF((float)x_avg, (float)Y));
|
}
|
else if(list.Count == 0)
|
{
|
//使用前一个点的值,X和Y都是前一个点的值
|
down_list.Add(new PointF((float)X, (float)Y));
|
up_list.Add(new PointF((float)X, (float)Y));
|
avg_list.Add(new PointF((float)X, (float)Y));
|
}
|
|
}
|
|
//CurveFitHelper h0 = new CurveFitHelper(down_list);
|
//var curve0 = h0.GetFitCurve(count);
|
//CurveFitHelper h1 = new CurveFitHelper(up_list);
|
//var curve1 = h1.GetFitCurve(count);
|
return new List<List<PointF>> { down_list, up_list, avg_list };
|
}
|
public static double[] GetConfidenceRange(List<double> data, double confidenceLevel = 90, bool mode = true)//true是中位数，置信区间是false
|
{
|
//如果model==true,这计算中位数,比例为confidenceLevel
|
if (mode)
|
{
|
|
data.Sort();
|
|
double lowerBound = 0;
|
double upperBound = 0;
|
//取百分位confidenceLevel/100和1-confidenceLevel/100.0的数
|
int upperIndex = (int)Math.Round(data.Count * confidenceLevel / 100.0);
|
int lowerIndex = (int)Math.Round(data.Count * (1 - confidenceLevel / 100.0));
|
if (upperIndex == data.Count)
|
{
|
upperIndex = data.Count - 1;
|
}
|
if (lowerIndex < 0)
|
{
|
lowerIndex = 0;
|
}
|
return new double[] { data[lowerIndex], data[upperIndex] };
|
|
|
}
|
|
else
|
{
|
// 使用循环计算数据的最小值和最大值
|
double minValue = double.MaxValue;
|
double maxValue = double.MinValue;
|
foreach (double value in data)
|
{
|
if (value < minValue)
|
{
|
minValue = value;
|
}
|
if (value > maxValue)
|
{
|
maxValue = value;
|
}
|
}
|
|
// 使用循环计算样本均值
|
double sum = 0;
|
foreach (double value in data)
|
{
|
sum += value;
|
}
|
double mean = sum / data.Count;
|
|
// 使用循环计算标准差
|
double sumOfSquares = 0;
|
foreach (double value in data)
|
{
|
sumOfSquares += Math.Pow(value - mean, 2);
|
}
|
double stdDev = Math.Sqrt(sumOfSquares / (data.Count - 1));
|
|
double zValue = 0;
|
double alpha = 1 - confidenceLevel / 100.0;
|
|
double zValueOneTail = GetZValueOneTail(alpha / 2);
|
double marginOfError = zValueOneTail * stdDev;
|
double lowerBound = mean - marginOfError;
|
double upperBound = mean + marginOfError;
|
|
return new double[] { lowerBound, upperBound };
|
}
|
|
|
|
}
|
static double[,] zTable = {
|
{0.0000, 0.0039, 0.0078, 0.0117, 0.0156, 0.0195, 0.0234, 0.0274, 0.0314, 0.0353},
|
{0.0392, 0.0432, 0.0471, 0.0510, 0.0550, 0.0589, 0.0628, 0.0668, 0.0708, 0.0749},
|
{0.0789, 0.0830, 0.0871, 0.0912, 0.0953, 0.0994, 0.1035, 0.1076, 0.1117, 0.1158},
|
{0.1190, 0.1231, 0.1272, 0.1314, 0.1357, 0.1398, 0.1441, 0.1483, 0.1525, 0.1568},
|
{0.1611, 0.1654, 0.1697, 0.1740, 0.1783, 0.1826, 0.1869, 0.1912, 0.1955, 0.1999},
|
{0.2042, 0.2085, 0.2129, 0.2172, 0.2216, 0.2260, 0.2304, 0.2348, 0.2392, 0.2436},
|
{0.2480, 0.2524, 0.2568, 0.2612, 0.2658, 0.2703, 0.2747, 0.2793, 0.2838, 0.2883},
|
{0.2929, 0.2974, 0.3020, 0.3066, 0.3112, 0.3159, 0.3205, 0.3251, 0.3299, 0.3346},
|
{0.3393, 0.3441, 0.3487, 0.3535, 0.3583, 0.3632, 0.3679, 0.3727, 0.3776, 0.3823},
|
{0.3872, 0.3920, 0.3968, 0.4017, 0.4066, 0.4115, 0.4164, 0.4214, 0.4264, 0.4313}
|
};
|
static double GetZValueOneTail(double alpha)
|
{
|
if (alpha == 0.5)
|
{
|
return 0;
|
}
|
else if (alpha > 0.5)
|
{
|
// Use symmetry property of the normal distribution
|
return -GetZValueOneTail(1 - alpha);
|
}
|
double[] alphaValues = { 0.1, 0.05, 0.025, 0.01, 0.005, 0.001 };
|
double[] zValues = { 1.28, 1.645, 1.96, 2.33, 2.575, 3.09 };
|
int n = alphaValues.Length;
|
int i = 0;
|
while (i < n && alpha < alphaValues[i])
|
{
|
i++;
|
}
|
|
if (i == 0)
|
{
|
return 0;
|
}
|
else if (i == n)
|
{
|
return zValues[n - 1];
|
}
|
|
double x0 = alphaValues[i - 1];
|
double x1 = alphaValues[i];
|
double y0 = zValues[i - 1];
|
double y1 = zValues[i];
|
|
double zValue = ((y1 - y0) / (x1 - x0)) * (alpha - x0) + y0;
|
|
return zValue;
|
}
|
|
//数据频率分布图的英文怎么说,
|
/// <summary>
|
/// 获取数据频率的分布图
|
/// </summary>
|
/// <param name="data"></param>
|
/// <param name="numIntervals"></param>
|
/// <param name="filterNum">以最小值作为过滤的系数</param>
|
/// <returns></returns>
|
public static List<DataPoint> GetFrequencyDistribution(List<Result> data, out List<DataPoint>[] datapoints, int numIntervals = 10,double filterNum=5 )
|
{
|
|
List<DataPoint> resultPoints = new List<DataPoint>();
|
datapoints= new List<DataPoint>[numIntervals];
|
for(int i=0;i<numIntervals;i++)
|
{
|
datapoints[i]= new List<DataPoint>();
|
}
|
// 设置置信水平
|
|
|
// 设置区间数量
|
|
|
////在data中过滤和最大的5%的数据
|
//data = data.Where(d => d <= data.Min() * filterNum).ToList();
|
|
|
|
// 计算数据的最小值和最大值
|
double minValue = data.Min(x=>x.ObjFunctionValue);
|
double maxValue = data.Max(x => x.ObjFunctionValue);
|
|
// 计算区间宽度
|
double intervalWidth = (maxValue - minValue) / numIntervals;
|
|
// 初始化频率分布数组
|
int[] freqDist = new int[numIntervals];
|
|
// 计算各个区间的频率
|
foreach (var xd in data)
|
{
|
var x = xd.ObjFunctionValue;
|
int index = (int)((x - minValue) / intervalWidth);
|
if (index < 0) // 数据小于最小值
|
{
|
freqDist[0]++;
|
datapoints[0].Add(new DataPoint() {XValue= x,Tag=xd });
|
}
|
else if (index >= numIntervals) // 数据大于最大值
|
{
|
freqDist[numIntervals - 1]++;
|
datapoints[numIntervals - 1].Add(new DataPoint() { XValue = x, Tag = xd });
|
}
|
else // 数据在区间内
|
{
|
freqDist[index]++;
|
datapoints[index].Add(new DataPoint() { XValue = x, Tag = xd });
|
}
|
}
|
double[] freqDistPersent =new double[numIntervals];
|
double sum = freqDist.Sum();
|
//计算freqDist的频率百分比
|
for (int i = 0; i < numIntervals; i++)
|
{
|
freqDistPersent[i] = (freqDist[i] / (double)sum );
|
}
|
|
//// 计算样本均值和标准差
|
//double mean = data.Average();
|
//double stdDev = Math.Sqrt(data.Select(x => (x - mean) * (x - mean)).Sum() / (data.Length - 1));
|
|
//double zValue = 0;
|
//double alpha = 1 - confidenceLevel / 100.0;
|
|
//double zValueOneTail = GetZValueOneTail(alpha / 2);
|
|
////zValue = mean + zValueOneTail * stdDev / Math.Sqrt(data.Length);
|
//// 计算置信区间的范围
|
////double z = GetZValue(confidenceLevel); //chart1.DataManipulator.Statistics.InverseNormalDistribution(1 - confidenceLevel / 200.0); // 根据置信水平查找标准正态分布的临界值
|
//double marginOfError = zValueOneTail * stdDev;/// Math.Sqrt(data.Length);
|
//double lowerBound = mean - marginOfError;
|
//double upperBound = mean + marginOfError;
|
for (int i = 0; i < numIntervals; i++)
|
{
|
double xValue = minValue + i * intervalWidth + intervalWidth / 2;
|
double yValue = freqDistPersent[i];
|
if (isValid(xValue) && isValid(yValue))
|
resultPoints.Add(new DataPoint(xValue, yValue));
|
}
|
|
|
return resultPoints;
|
|
|
|
}
|
|
static bool isValid(double value)
|
{
|
if (value == (int)value && value < 0) return false;
|
return !Double.IsNaN(value) && !Double.IsInfinity(value);
|
}
|
|
|
}
|
|
//public static class ConfidenceRangeHelper
|
//{
|
|
|
//}
|
}
|
