Solved – Implementation of M-spline in R

I am implementing M-spline in R as defined here : http://www.fon.hum.uva.nl/praat/manual/spline.html and originally in Ramsay (1988).

In short, we define a list of knots \$t\$ such that :

\$0=t_1=…=t_k<t_{k+1}<…<t_{k+m}<t_{k+m+1}=…=t_{k+m+k}=1\$

So the \$k\$ first knots are 0 and the \$k\$ last knots are 1 with \$m\$ inner knots.

The Mi spline of order k with knots t is define recursively as below :

\$M_i(x|1,t) = 1 / (t_{i+1} – t_i), t_i ≤ x < t_{i+1},\$ 0 otherwise

\$M_i(x|k,t) = k [(x–t_i)M_i(x|k–1,t) + (t_{i+k}–x)M_{i+1}(x|k–1,t)] / ((k–1)(t_{i+k}–t_i))\$

The problem is that if I want the M-splines of order 3 (k=3), the \$M_1\$ spline and the \$M_{k+m}\$ spline do not exist because when it calculate \$M_i(x|k–1,t)\$ it divides by \$(t_{i+k}–t_i)\$ but since now \$k=3-1=2\$ and \$i=1\$ (or \$m+k\$), \$(t_{i+k}–t_i)\$ is \$(t_{1+2}–t_1)=0\$ since the first three knots are 0 (same thing for the last one).

So the definition seems to have a problem. Here is my implementation in R with an example of the problem. Here k=3 and m=3.

``ts = c(0,0,0,0.3,0.5,0.6,1,1,1) Mk = function(i,k,x,ts){   if(k==1){     if(ts[i]<=x && x<ts[i+1]){1/(ts[i+1]-ts[i])}     else{0}   }   else{     #print(paste(i,k))     #print((ts[i+k]-ts[i]))     k*((x-ts[i])*Mk(i,k-1,x,ts)+(ts[i+k]-x)*Mk(i+1,k-1,x,ts))/((k-1)*(ts[i+k]-ts[i]))    } } Mk(1,3,.4,ts) ``

Note that my code and the definition seems to work for all the "inner" splines.

Thanks!

Contents

``    if(ts[i+k]-ts[i]==0){0} ``