How to obrain p-values for coefficients with polyval - polyparci?

I am using polyfit / polyval to obtain coefficients for a simple quadratic relationship:

[p,S] = polyfit(x2,y2,order);

[y_fit,delta] = polyval(p,x2,S);

With polyparci i can compute the confidence intervals fo the coefficients:

ci = polyparci(p,S);

but there is no built-in way to have p-values for those coefficient? I know it's implicit in polyparci but I can't figure out a simple way of extracting that specific output.

Thanks in advance

[Torsten]

What do you mean by "obtain p values" ? You prescribe the p values and MATLAB returns the confidence intervals according to your p-value. Or do you want to prescribe confidence intervals to get the p values ?

[the cyclist]

@Torsten, I think you may be confusing p-value with what is usually called "alpha", the Type 1 error rate threshold. alpha is prescribed ahead of modeling, and will be a factor in calculating the confidence interval. p-value is an output of the model.

The null hypothesis can be rejected if the p-value is smaller than the prescribed alpha.

[Star Strider]

The ‘polyparci’ function returns the 95% parameter confidence limits. or ‘alpha’ if it is provided. The parameters are ‘significant’ at ‘p=0.05’ (or whatever ‘alpha’ value is chosen) if they have the same signs. If the signs differ, then the confidence intervals include 0, and that parameter is likely not necessary in the regression.

Also, I just made a major upgrade to ‘polyparci’. See my Answer for details.

EDIT — (10 Jan 2024 at 18:30)

Answer 1

[Sulaymon Eshkabilov]

Note that polyparci() is not a MATLAB built in fcn and it is posted here. Note that it is quite stratforward to obtain p-values of a polynomial using fitlm(). Here is a short demo example for a cubic polynomial fit:

x = -2:.25:2;
y = 3*x.^2-5*x+10+randn(size(x))*10;
%% fitlm() to compute p-values of fit model coefficients:
fitlm(x,y, 'poly3')

ans = 
Linear regression model:
    y ~ 1 + x1 + x1^2 + x1^3

Estimated Coefficients:
                   Estimate      SE       tStat       pValue  
                   ________    ______    _______    __________

    (Intercept)     10.133     2.5776     3.9313     0.0017214
    x1             -17.456     3.5409    -4.9299    0.00027507
    x1^2            4.4398     1.2838     3.4584     0.0042396
    x1^3            4.0843     1.2104     3.3744     0.0049814


Number of observations: 17, Error degrees of freedom: 13
Root Mean Squared Error: 7.06
R-squared: 0.775,  Adjusted R-Squared: 0.723
F-statistic vs. constant model: 14.9, p-value = 0.000167

% Here is how to display confidence interval of 95%

[FM, S]=polyfit(x,y,3);                 % Fit Model coeffs: FM
[FM_vals, delta] = polyval(FM,x, S);    % Prediction Interval at 95%: delta
plot(x,y, 'ko', 'MarkerFaceColor','c', 'MarkerSize', 7)
hold on
plot(x, FM_vals, 'k-', 'linewidth',2)
plot(x, FM_vals+2*delta, 'r--', x, FM_vals-2*delta, 'r--')
legend('Data', 'Fit Model', '95% Prediction Interval')
grid on
xlabel('x')
ylabel('y(x)')

[the cyclist]

@LEONARDO CIAMBEZI, note that fitlm is a function from the Statistics Learning Toolbox, not base MATLAB.

@Sulaymon Eshkabilov, note that your fitlm example is using 3rd-order polynomial to fit a 2nd-order function, so your fitting coefficients are somewhat nonsensical. A more suitable example is

x = -2:.25:2;
y = 3*x.^2-5*x+10+randn(size(x))*0.4;
%% fitlm() to compute p-values of fit model coefficients:
fitlm(x,y, 'poly2')

ans = 
Linear regression model:
    y ~ 1 + x1 + x1^2

Estimated Coefficients:
                   Estimate       SE        tStat       pValue  
                   ________    ________    _______    __________

    (Intercept)     10.132      0.12698     79.797     5.128e-20
    x1             -5.0373     0.068916    -73.093    1.7471e-19
    x1^2             2.913     0.063242     46.062    1.0926e-16


Number of observations: 17, Error degrees of freedom: 14
Root Mean Squared Error: 0.348
R-squared: 0.998,  Adjusted R-Squared: 0.998
F-statistic vs. constant model: 3.73e+03, p-value = 8.06e-20

[Chinmayee L M]

Can I use fitlm to fit the equation 'y = a*exp(b*x)' ? If yes, how?

The fit function from the curve fit toolbox fits this equation with 'exp1' model name. But, I don't get p-values from it.

[Star Strider]

‘Can I use fitlm to fit the equation 'y = a*exp(b*x)' ? If yes, how?’

No. That requires fitnlm.

Example —

x = (0:0.1:20).';
y = 5*exp(-0.25*x) + randn(size(x));

mdl = fitnlm(x, y, @(b,x) b(1).*exp(b(2).*x), randn(2,1))

mdl = 
Nonlinear regression model:
    y ~ b1*exp(b2*x)

Estimated Coefficients:
          Estimate      SE        tStat       pValue  
          ________    _______    _______    __________

    b1     5.2145     0.29076     17.934    1.9789e-43
    b2     -0.263     0.02103    -12.506    7.3825e-27


Number of observations: 201, Error degrees of freedom: 199
Root Mean Squared Error: 0.92
R-Squared: 0.672,  Adjusted R-Squared 0.67
F-statistic vs. zero model: 314, p-value = 2.96e-62

[ypred,yci] = predict(mdl, x);

figure
hp1 = plot(x, y, 'p', 'DisplayName','Data');
hold on
hp2 = plot(x, ypred, '-r', 'DisplayName','Regression');
hp3 =  plot(x, yci, '--r', 'DisplayName', '95% Confidence Interval');
hold off
grid
legend([hp1; hp2; hp3(1)], 'Location','best')

.

[the cyclist]

The equation y = a*exp(b*x) is non-linear in the parameter b, so you should not try to fit a linear model to it.

But, if you take the log of both sides of the equation, you get log(y) = log(a) + b*x, which you can sensibly fit with a linear model.

Using an analogous example to @Star Strider's,

x = (0:0.1:5).';
y = 5*exp(-0.25*x) + 0.1*randn(size(x));

logy = log(y); % Transformed data

% Non-linear model
mdl_nlm    = fitnlm(x, y, @(b,x) b(1).*exp(b(2).*x), randn(2,1))

mdl_nlm = 
Nonlinear regression model:
    y ~ b1*exp(b2*x)

Estimated Coefficients:
          Estimate       SE        tStat       pValue  
          ________    ________    _______    __________

    b1      4.9925      0.0387        129    1.0329e-63
    b2    -0.24857    0.003908    -63.606    9.2785e-49


Number of observations: 51, Error degrees of freedom: 49
Root Mean Squared Error: 0.108
R-Squared: 0.99,  Adjusted R-Squared 0.989
F-statistic vs. zero model: 2.03e+04, p-value = 3.17e-72

% Linear model on log-transformed equation
mdl_log_lm = fitlm(x, logy)

mdl_log_lm = 
Linear regression model:
    y ~ 1 + x1

Estimated Coefficients:
                   Estimate       SE         tStat       pValue  
                   ________    _________    _______    __________

    (Intercept)      1.6089     0.013667     117.72    9.0257e-62
    x1             -0.24948    0.0047108    -52.959    6.4673e-45


Number of observations: 51, Error degrees of freedom: 49
Root Mean Squared Error: 0.0495
R-squared: 0.983,  Adjusted R-Squared: 0.982
F-statistic vs. constant model: 2.8e+03, p-value = 6.47e-45

Note that Star Strider's and my second coefficient are equal (to within the model tolerance), and my first coefficient is equal to the log of his (again with a bit of tolerance):

log(5.045)

ans = 1.6184

Answer 2

[Star Strider]

‘... but there is no built-in way to have p-values for those coefficient?’

There is now. I did not realise that there was a need for it, however as an improvement to polyparci, I added that as the second output (along with some significant improvements to it) and uploaded it about a week ago.