10
$\begingroup$

I need to fit a spline function to a data set. I tried with bs, ns and smooth.spline. In my case the curve obtained by smooth.spline follows the trend in data better than with the bs and ns.

However, I do not know how to

  1. obtain the equation of the fitted smooth.spline function?
    (namely, how to interpret the equation coefficients of the smooth.spline output)
  2. compute the $R^2$ of the fitted function?
  3. calculate the analytical first derivative of the fitted function, and obtain its equation (with actual coefficients)?

An example of data set can be computed as

library(stats)
x <- 1:11
y <- c(0.2,0.40, 0.6, 0.75, 0.88, 0.99, 1.1, 1.15, 1.16, 1.16, 1.16 )
plot(x, y)
spline <- smooth.spline(x,y)
lines(spline)
$\endgroup$
6
$\begingroup$

This fits a natural spline (linear tail restricted) using the truncated power basis. In this example, default knots (based on quantiles of the predictor) are not used; instead we specify 4 knots. The only way to get a test of goodness of fit is to postulate a richer model than this and see if it improves the model fitted below. But anova() tests the goodness of fit of a linear relationship by pooling the nonlinear terms into a composite ("chunk") test (F=175.38).

require(rms)
x <- 1:11
y <- c(0.2,0.40, 0.6, 0.75, 0.88, 0.99, 1.1, 1.15, 1.16, 1.16, 1.16 )
dd <- datadist(x); options(datadist='dd')

f <- ols(y ~ rcs(x, c(3, 5, 7, 9)))
f

Linear Regression Model

ols(formula = y ~ rcs(x, c(3, 5, 7, 9)))

                Model Likelihood     Discrimination    
                   Ratio Test           Indexes        
Obs       11    LR chi2     66.08    R2       0.998    
sigma 0.0201    d.f.            3    R2 adj   0.996    
d.f.       7    Pr(> chi2) 0.0000    g        0.383    

Residuals

      Min        1Q    Median        3Q       Max 
-0.027360 -0.011739  0.001227  0.009892  0.031166 

          Coef    S.E.   t     Pr(>|t|)
Intercept  0.0465 0.0224  2.08 0.0762  
x          0.1741 0.0072 24.18 <0.0001 
x'        -0.1004 0.0311 -3.23 0.0144  
x''        0.0542 0.0913  0.59 0.5715  

anova(f)

                Analysis of Variance          Response: y 

 Factor     d.f. Partial SS  MS           F      P     
 x          3    1.152321844 0.3841072814 946.15 <.0001
  Nonlinear 2    0.142398208 0.0711991039 175.38 <.0001
 REGRESSION 3    1.152321844 0.3841072814 946.15 <.0001
 ERROR      7    0.002841792 0.0004059703              

ggplot(Predict(f)) + geom_point(aes(x=x, y=y), data=data.frame(x,y))
Function(f)   ## if have latex installed can also use latex(f)

function(x = 6) {0.046475489+0.17411942* x-0.002790266*pmax(x-    3,0)^3+0.0015048699*pmax(x-5,0)^3+0.0053610582*pmax(x-7,0)^3-0.0040756621*pmax(x-9,0)^3 }

Function re-expresses the restricted cubic spline in simplest form. The first derivative is:

function(x) 0.174 - 3 * 0.00279 * pmax(x - 3, 0) ^ 2 + 3 * 0.0015 * pmax(x - 5, 0) ^ 2 + ...
| cite | improve this answer | |
$\endgroup$
0
$\begingroup$

I would suggest using interpSpline from the splines package. Building on your example, you could use:

library(splines)
x <- 1:11
y <- c(0.2,0.40, 0.6, 0.75, 0.88, 0.99, 1.1, 1.15, 1.16, 1.16, 1.16 )
spline <- interpSpline(x,y)
plot(spline)
points(x,y)

The object spline then contains the coefficients you wanted as a part of your 1. and 3. question, getting the equations of the first derivative by differentiating the original equations constructed using these coefficients.

As for your second question, I am not sure what you are trying to do. If you are simply interpolating a cubic spline through all of your data points, i.e. using all of them as knots, as you did in your example, then the curve fits the points perfectly.

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ Thank very much you Fato39 for your answer! However I'd like to specify a fix small number of knots (e.g. 3) belonging to the spline, which is possible with smooth.spline. Additionally I'd also like to set at which x values these knots occur. Among the different spline functions present in R, I have not found one where I can: i) fit the spline while deciding how many knots and where along the x-axis; ii) obtaining the equation of the fitted spline; iii) obtain the derivative of the fitted function and iv) possibly having a statistics the goodness of fit of the spline. $\endgroup$ – Francesco Apr 27 '15 at 17:09

Not the answer you're looking for? Browse other questions tagged or ask your own question.