Add comprehensive statistical metrics to NonlinearMinimizationResult

- TStatistics and PValues for parameter significance testing
- ConfidenceIntervalHalfWidths for parameter precision
- Dependencies to measure parameter correlations
- Goodness-of-fit statistics

Author: diluculo

src/Numerics.Tests/OptimizationTests/NonLinearCurveFittingTests.cs

@@ -660,6 +660,36 @@ public void PollutionWithWeights()
             {
                 AssertHelpers.AlmostEqualRelative(PollutionBest[i], result.MinimizingPoint[i], 4);
             }
+
+            // Check statistics
+            AssertHelpers.AlmostEqualRelative(0.908, result.RSquared, 2);
+            AssertHelpers.AlmostEqualRelative(0.885, result.AdjustedRSquared, 2);
+            AssertHelpers.AlmostEqualRelative(24.0096, result.StandardError, 2);
+
+            // Check parameter statistics (using expected values)
+            var expectedStdErrors = new double[] { 10.7, 0.064296 };
+            var expectedTStats = new double[] { 21.045, 6.2333 };
+            var expectedPValues = new double[] { 3.0134e-05, 0.0033745 };
+            // Expected confidence interval bounds
+            var expectedCI = new double[,]
+                {
+                    { 195.4650, 254.8788 },
+                    { 0.2223, 0.5793 }
+                };
+
+            for (var i = 0; i < result.MinimizingPoint.Count; i++)
+            {
+                AssertHelpers.AlmostEqualRelative(expectedStdErrors[i], result.StandardErrors[i], 3);
+                AssertHelpers.AlmostEqualRelative(expectedTStats[i], result.TStatistics[i], 3);
+                AssertHelpers.AlmostEqualRelative(expectedPValues[i], result.PValues[i], 3);
+
+                // Calculate and check confidence interval bounds
+                var lowerBound = result.MinimizingPoint[i] - result.ConfidenceIntervalHalfWidths[i];
+                var upperBound = result.MinimizingPoint[i] + result.ConfidenceIntervalHalfWidths[i];
+
+                AssertHelpers.AlmostEqualRelative(expectedCI[i, 0], lowerBound, 3);
+                AssertHelpers.AlmostEqualRelative(expectedCI[i, 1], upperBound, 3);
+            }
         }
 
         #endregion Weighted Nonlinear Regression

src/Numerics/Optimization/NonlinearMinimizationResult.cs

@@ -1,9 +1,19 @@
-using MathNet.Numerics.LinearAlgebra;
+using MathNet.Numerics.Distributions;
+using MathNet.Numerics.LinearAlgebra;
+using System;
+using System.Linq;
 
 namespace MathNet.Numerics.Optimization
 {
+    /// <summary>
+    /// Represents the result of a nonlinear minimization operation, including
+    /// the optimal parameters and various statistical measures of fitness.
+    /// </summary>
     public class NonlinearMinimizationResult
     {
+        /// <summary>
+        /// The objective model evaluated at the minimum point.
+        /// </summary>
         public IObjectiveModel ModelInfoAtMinimum { get; }
 
         /// <summary>
@@ -16,6 +26,30 @@ public class NonlinearMinimizationResult
         /// </summary>
         public Vector<double> StandardErrors { get; private set; }
 
+        /// <summary>
+        /// Returns the t-statistics for each parameter (parameter value / standard error).
+        /// These measure how many standard deviations each parameter is from zero.
+        /// </summary>
+        public Vector<double> TStatistics { get; private set; }
+
+        /// <summary>
+        /// Returns the p-values for each parameter based on t-distribution.
+        /// Lower p-values indicate statistically significant parameters.
+        /// </summary>
+        public Vector<double> PValues { get; private set; }
+
+        /// <summary>
+        /// Returns the dependency values for each parameter, measuring how linearly related
+        /// each parameter is to the others. Values close to 1 indicate high dependency (multicollinearity).
+        /// </summary>
+        public Vector<double> Dependencies { get; private set; }
+
+        /// <summary>
+        /// Returns the half-width of the confidence intervals for each parameter at the 95% confidence level.
+        /// These represent the margin of error for each parameter estimate.
+        /// </summary>
+        public Vector<double> ConfidenceIntervalHalfWidths { get; private set; }
+
         /// <summary>
         /// Returns the y-values of the fitted model that correspond to the independent values.
         /// </summary>
@@ -31,19 +65,69 @@ public class NonlinearMinimizationResult
         /// </summary>
         public Matrix<double> Correlation { get; private set; }
 
+        /// <summary>
+        /// Returns the residual sum of squares at the minimum point, 
+        /// calculated as F(p) = 1/2 * ∑(residuals²).
+        /// </summary>
+        public double MinimizedValue => ModelInfoAtMinimum.Value;
+
+        /// <summary>
+        /// Root Mean Squared Error (RMSE) - measures the average magnitude of errors.
+        /// </summary>
+        public double RootMeanSquaredError { get; private set; }
+
+        /// <summary>
+        /// Coefficient of determination (R-Square) - proportion of variance explained by the model.
+        /// </summary>
+        public double RSquared { get; private set; }
+
+        /// <summary>
+        /// Adjusted R-Squre - accounts for the number of predictors in the model.
+        /// </summary>
+        public double AdjustedRSquared { get; private set; }
+
+        /// <summary>
+        /// Residual standard error of the regression, calculated as sqrt(RSS/df) where RSS is the 
+        /// residual sum of squares and df is the degrees of freedom. This measures the average 
+        /// distance between the observed values and the fitted model.
+        /// </summary>
+        public double StandardError { get; private set; }
+
+        /// <summary>
+        /// Pearson correlation coefficient between observed and predicted values.
+        /// </summary>
+        public double CorrelationCoefficient { get; private set; }
+
+        /// <summary>
+        /// Number of iterations performed during optimization.
+        /// </summary>
         public int Iterations { get; }
 
+        /// <summary>
+        /// Reason why the optimization algorithm terminated.
+        /// </summary>
         public ExitCondition ReasonForExit { get; }
 
+        /// <summary>
+        /// Creates a new instance of the NonlinearMinimizationResult class.
+        /// </summary>
+        /// <param name="modelInfo">The objective model at the minimizing point.</param>
+        /// <param name="iterations">The number of iterations performed.</param>
+        /// <param name="reasonForExit">The reason why the algorithm terminated.</param>
         public NonlinearMinimizationResult(IObjectiveModel modelInfo, int iterations, ExitCondition reasonForExit)
         {
             ModelInfoAtMinimum = modelInfo;
             Iterations = iterations;
             ReasonForExit = reasonForExit;
 
             EvaluateCovariance(modelInfo);
+            EvaluateGoodnessOfFit(modelInfo);
         }
 
+        /// <summary>
+        /// Evaluates the covariance matrix, correlation matrix, standard errors, t-statistics and p-values.
+        /// </summary>
+        /// <param name="objective">The objective model at the minimizing point.</param>
         void EvaluateCovariance(IObjectiveModel objective)
         {
             objective.EvaluateAt(objective.Point); // Hessian may be not yet updated.
@@ -66,15 +150,208 @@ void EvaluateCovariance(IObjectiveModel objective)
             {
                 StandardErrors = Covariance.Diagonal().PointwiseSqrt();
 
+                // Calculate t-statistics and p-values
+                TStatistics = Vector<double>.Build.Dense(StandardErrors.Count);
+                for (var i = 0; i < StandardErrors.Count; i++)
+                {
+                    TStatistics[i] = StandardErrors[i] > double.Epsilon
+                        ? objective.Point[i] / StandardErrors[i]
+                        : double.NaN;
+                }
+
+                // Calculate p-values based on t-distribution with DegreeOfFreedom
+                PValues = Vector<double>.Build.Dense(TStatistics.Count);
+                var tDist = new StudentT(0, 1, objective.DegreeOfFreedom);
+
+                // Calculate confidence interval half-widths (for 95% confidence)
+                ConfidenceIntervalHalfWidths = Vector<double>.Build.Dense(StandardErrors.Count);
+                var tCritical = tDist.InverseCumulativeDistribution(0.975); // Two-tailed, 95% confidence
+
+                for (var i = 0; i < TStatistics.Count; i++)
+                {
+                    // Two-tailed p-value from t-distribution
+                    var tStat = Math.Abs(TStatistics[i]);
+                    PValues[i] = 2 * (1 - tDist.CumulativeDistribution(tStat));
+
+                    // Calculate confidence interval half-width
+                    ConfidenceIntervalHalfWidths[i] = tCritical * StandardErrors[i];
+                }
+
                 var correlation = Covariance.Clone();
                 var d = correlation.Diagonal().PointwiseSqrt();
                 var dd = d.OuterProduct(d);
                 Correlation = correlation.PointwiseDivide(dd);
+
+                // Calculate dependencies (measure of multicollinearity)
+                Dependencies = Vector<double>.Build.Dense(Correlation.RowCount);
+                for (var i = 0; i < Correlation.RowCount; i++)
+                {
+                    // For each parameter, calculate 1 - 1/VIF where VIF is the variance inflation factor
+                    // VIF is calculated as 1/(1-R²) where R² is the coefficient of determination when
+                    // parameter i is regressed against all other parameters
+
+                    // Extract the correlation coefficients related to parameter i (excluding self-correlation)
+                    var correlations = Vector<double>.Build.Dense(Correlation.RowCount - 1);
+                    var index = 0;
+                    for (var j = 0; j < Correlation.RowCount; j++)
+                    {
+                        if (j != i)
+                        {
+                            correlations[index++] = Correlation[i, j];
+                        }
+                    }
+
+                    // Calculate the square of the multiple correlation coefficient
+                    // In the simple case, this is the maximum of the squared individual correlations
+                    var maxSquaredCorrelation = 0.0;
+                    for (var j = 0; j < correlations.Count; j++)
+                    {
+                        var squaredCorr = correlations[j] * correlations[j];
+                        if (squaredCorr > maxSquaredCorrelation)
+                        {
+                            maxSquaredCorrelation = squaredCorr;
+                        }
+                    }
+
+                    // Calculate dependency = 1 - 1/VIF
+                    var maxSquaredCorrelationCapped = Math.Min(maxSquaredCorrelation, 0.9999);
+                    var vif = 1.0 / (1.0 - maxSquaredCorrelationCapped);
+                    Dependencies[i] = 1.0 - 1.0 / vif;
+                }
             }
             else
             {
                 StandardErrors = null;
+                TStatistics = null;
+                PValues = null;
+                ConfidenceIntervalHalfWidths = null;
                 Correlation = null;
+                Dependencies = null;
+            }
+        }
+
+        /// <summary>
+        /// Evaluates goodness of fit statistics like R-squared, RMSE, etc.
+        /// </summary>
+        /// <param name="objective">The objective model at the minimizing point.</param>
+        void EvaluateGoodnessOfFit(IObjectiveModel objective)
+        {
+            // Note: GoodnessOfFit class methods do not support weights, so we apply weighting manually here.
+
+            // Check if we have the essentials for calculating statistics
+            var hasResiduals = objective.Residuals != null;
+            var hasObservations = objective.ObservedY != null;
+            var hasModelValues = objective.ModelValues != null;
+            var hasSufficientDof = objective.DegreeOfFreedom >= 1;
+
+            // Set values to NaN if we can't calculate them
+            RootMeanSquaredError = double.NaN;
+            RSquared = double.NaN;
+            AdjustedRSquared = double.NaN;
+            StandardError = double.NaN;
+            CorrelationCoefficient = double.NaN;
+
+            // Need residuals and sufficient DOF for most calculations
+            if (!hasResiduals || !hasSufficientDof)
+            {
+                return;
+            }
+
+            var n = hasObservations ? objective.ObservedY.Count : objective.Residuals.Count;
+            var dof = objective.DegreeOfFreedom;
+
+            // Calculate sum of squared residuals using vector operations
+            var ssRes = objective.Residuals.DotProduct(objective.Residuals);
+
+            // Guard against zero or negative SSR
+            if (ssRes <= 0)
+            {
+                RootMeanSquaredError = 0;
+                StandardError = 0;
+
+                // Only calculate these if we have observations
+                if (hasObservations)
+                {
+                    RSquared = 1.0;
+                    AdjustedRSquared = 1.0;
+                }
+
+                // Only calculate if we have model values and observations
+                if (hasModelValues && hasObservations)
+                {
+                    CorrelationCoefficient = 1.0;
+                }
+                return;
+            }
+
+            // Calculate standard error and RMSE, which only require residuals
+            StandardError = Math.Sqrt(ssRes / dof);
+            RootMeanSquaredError = Math.Sqrt(ssRes / n);
+
+            // The following statistics require observations
+            if (!hasObservations)
+            {
+                return;
+            }
+
+            // Calculate total sum of squares
+            double ssTot;
+
+            // If weights are present, calculate weighted total sum of squares
+            if (objective.Weights != null)
+            {
+                var weightSum = 0.0;
+                var weightedSum = 0.0;
+
+                for (var i = 0; i < n; i++)
+                {
+                    var weight = objective.Weights[i, i];
+                    weightSum += weight;
+                    weightedSum += weight * objective.ObservedY[i];
+                }
+
+                // Avoid division by zero
+                var weightedMean = weightSum > double.Epsilon ? weightedSum / weightSum : 0;
+
+                ssTot = 0.0;
+                for (var i = 0; i < n; i++)
+                {
+                    var weight = objective.Weights[i, i];
+                    var dev = objective.ObservedY[i] - weightedMean;
+                    ssTot += weight * dev * dev;
+                }
+            }
+            else
+            {
+                // Unweighted case - use vector operations for total sum of squares
+                var yMean = objective.ObservedY.Average();
+                var deviations = objective.ObservedY.Subtract(yMean);
+                ssTot = deviations.DotProduct(deviations);
+            }
+
+            // Guard against zero or negative total sum of squares
+            if (ssTot <= double.Epsilon)
+            {
+                RSquared = 0.0;
+                AdjustedRSquared = 0.0;
+            }
+            else
+            {
+                // Calculate R-squared directly
+                RSquared = 1 - (ssRes / ssTot);
+
+                // Calculate adjusted R-squared using the ratio of mean squares
+                AdjustedRSquared = 1 - (ssRes / dof) / (ssTot / (n - 1));
+
+                // Ensure adjusted R-squared is not greater than 1 or less than 0
+                AdjustedRSquared = Math.Min(1.0, Math.Max(0.0, AdjustedRSquared));
+            }
+
+            // Only calculate correlation coefficient if we have model values
+            if (hasModelValues)
+            {
+                // Calculate correlation coefficient between observed and predicted
+                CorrelationCoefficient = GoodnessOfFit.R(objective.ModelValues, objective.ObservedY);
             }
         }
     }