[R] 1st derivative {mgcv} gam smooth

[Henric Nilsson]

On On, 2005-11-23, 05:43, Tomas skrev:
> Dear R-hep,
>           I'm trying to get the first derivative of a smooth from a gam
> model like:
> model<-gam(y~s(x,bs="cr", k=5)+z) and need the derivative: ds(x)/dx. Since
> coef(model) give me all the parameters, including the parameters of the
> basis, I just need the 1st derivative of the basis s(x).1, s(x).2, s(x).3,
> s(x).4. If the basis were generated with the function spline.des() from
> the
> spline library I would be able to get the derivative of the basis and
> would
> be pretty much done. (spline.des() has the option of computing 1st
> derivative)
>
> The difficulty I'm having is that I can't get how the basis were
> generated.
> If I define an "s" object like: object<-s(x,bs="cr",k=5) and use
> smooth.construct() to generate the basis, it gives 5 basis but the gam
> object "model" gives only 4 s(x).i elements. It generalizes for different

Try

> model2 <- gam(y ~ s(x, bs = "cr", k = 5) + z, control =
gam.control(absorb.cons = FALSE))
> coef(model2)
(Intercept)         z     s(x).1     s(x).2     s(x).3    s(x).4    s(x).5
  43.675698 -2.685739 -43.305376 -34.325047 -14.015483 20.155349 71.490557

See ?gam.control

HTH,
Henric


> values of k. See the example:
>
> library(mgcv)
> y<-c(1,5,6,12,13,14,45,42,56,68,89,120)
> z<-c(rep(0,5),rep(1,7))
> x<-c(seq(1:12))
> data<-data.frame(y,x,z)
>
> model<-gam(y~s(x,bs="cr", k=5)+z)
>
> object<-s(x,bs="cr", k=5)
> basis<-smooth.construct(object,data,knots=NULL)
>
> summary(model)
> .
> .
> .
> Parametric coefficients:
>             Estimate Std. Error t value Pr(>|t|)
> (Intercept)   43.676      5.492   7.953 6.66e-05 ***
> z             -2.686      9.041  -0.297    0.775
> ---
> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>
> Approximate significance of smooth terms:
>       edf Est.rank     F  p-value
> s(x) 2.53     4.00 35.18 6.48e-05 ***
>
> coef(model)
> (Intercept)           z      s(x).1      s(x).2      s(x).3      s(x).4
>  43.6756983  -2.6857394 -20.9429503  -0.6333857  33.5374463  84.8726537
>
> dim(basis$X)
> [1] 12  5
>
> #Here are the basis
> basis$X
>
>               [,1]        [,2]        [,3]        [,4]         [,5]
>  [1,]  1.000000000  0.00000000  0.00000000  0.00000000  0.000000000
>  [2,]  0.551840721  0.55522164 -0.13523666  0.03380917 -0.005634861
>  [3,]  0.180959536  0.93527960 -0.14682838  0.03670709 -0.006117849
>  [4,] -0.035821616  0.96624450  0.08768917 -0.02173446  0.003622411
>  [5,] -0.076875604  0.62353762  0.54792315 -0.11350220  0.018917033
>  [6,] -0.027771815  0.17264141  0.93603091 -0.09708061  0.016180101
>  [7,]  0.016180101 -0.09708061  0.93603091  0.17264141 -0.027771815
>  [8,]  0.018917033 -0.11350220  0.54792315  0.62353762 -0.076875604
>  [9,]  0.003622411 -0.02173446  0.08768917  0.96624450 -0.035821616
> [10,] -0.006117849  0.03670709 -0.14682838  0.93527960  0.180959536
> [11,] -0.005634861  0.03380917 -0.13523666  0.55522164  0.551840721
> [12,]  0.000000000  0.00000000  0.00000000  0.00000000  1.000000000
>
> #and their knots
>
> basis$xp
> [1]  1.00  3.75  6.50  9.25 12.00
>
> 	[[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide!
> http://www.R-project.org/posting-guide.html