gramian: Gramian matrix, norms, and inner products of splines
In Splinets: Functional Data Analysis using Splines and Orthogonal Spline Bases

gramian

R Documentation

Gramian matrix, norms, and inner products of splines

Description

The function performs evaluation of the matrix of the inner products \int S(t) * T(t) dt of all the pairs of splines S, T from the input object. The program utilizes the Taylor expansion of splines, see the reference for details.

Usage

gramian(Sp, norm_only = FALSE, sID = NULL, Sp2 = NULL, s2ID = NULL)

Arguments

`Sp`	`Splinets` object;
`norm_only`	logical, indicates if only the square norm of the elements in the input object is calculated; The default is `norm_only=FALSE`;
`sID`	vector of integers, the indicies specifying splines in the `Splinets` list `Sp` to be evaluated; If `sID=NULL` (default), then the inner products for all the pairs taken from the object are evaluated.
`Sp2`	`Splinets` object, the optional second `Splinets`-object; The inner products between splines in `Sp` and in `Sp2` are evaluated, i.e. the cross-gramian matrix.
`s2ID`	vector of integers, the indicies specifying splines in the `Sp2` to be considered in the cross-gramian;

Details

If there is only one input Splinet-object, then the non-negative symmetrix matrix of the splines in this object is returned. If there are two input Splinet-objects, then the m x r matrix of the cross-inner product is returned, where m is the number of splines in the first object and r is their number in the second one. If only the norms are evaluated (norm_only= TRUE) it is always evaluating the norms of the first object. In the case of two input Splinets-objects, they should be over the same set of knots and of the same smoothness order.

Value

norm_only=FALSE – the Gram matrix of inner products of the splines within the input Splinets-objects is returned,
Sp2 = NULL – the non-negative definite matrix of the inner products of splines in Sp is returned,
both Sp and Sp2 are non-NULL and contain splines S_i's and T_j's, respectively == the cross-gramian matris of the inner products for the pairs of splines (S_i,T_j) is returned,
norm_only=FALSE– the vector of the norms of Sp is returned.

References

Liu, X., Nassar, H., Podgorski, K. "Dyadic diagonalization of positive definite band matrices and efficient B-spline orthogonalization." Journal of Computational and Applied Mathematics (2022) <https://doi.org/10.1016/j.cam.2022.114444>.

Podgorski, K. (2021) "Splinets – splines through the Taylor expansion, their support sets and orthogonal bases." <arXiv:2102.00733>.

Nassar, H., Podgorski, K. (2023) "Splinets 1.5.0 – Periodic Splinets." <arXiv:2302.07552>

Examples

#---------------------------------------#
#---- Simple three splines example -----# 
#---------------------------------------#
n=25; k=3
xi=sort(runif(n+2)); xi[1]=0; xi[n+2]=1
#Defining support ranges for three splines
supp=matrix(c(2,12,4,20,6,25),byrow=TRUE,ncol=2)
#Initial random matrices of the derivative for each spline
SS1=matrix(rnorm((supp[1,2]-supp[1,1]+1)*(k+1)),ncol=(k+1)) 
SS2=matrix(rnorm((supp[2,2]-supp[2,1]+1)*(k+1)),ncol=(k+1)) 
SS3=matrix(rnorm((supp[3,2]-supp[3,1]+1)*(k+1)),ncol=(k+1)) 
spl=construct(xi,k,SS1,supp[1,]) #constructing the first correct spline
nspl=construct(xi,k,SS2,supp[2,])
spl=gather(spl,nspl) #the second and the first ones
nspl=construct(xi,k,SS3,supp[3,])
spl=gather(spl,nspl) #the third is added

plot(spl)
gramian(spl)
gramian(spl, norm_only = TRUE)
gramian(spl, sID = c(1,3))
gramian(spl,sID=c(2,3),Sp2=spl,s2ID=c(1)) #the cross-Gramian matrix  

#-----------------------------------------#
#--- Example with varying support sets ---#
#-----------------------------------------#
n=40; xi=seq(0,1,by=1/(n+1)); k=2; 
support=list(matrix(c(2,9,15,24,30,37),ncol=2,byrow = TRUE))
sp=new("Splinets",knots=xi,smorder=k,supp=support) 
m=sum(sp@supp[[1]][,2]-sp@supp[[1]][,1]+1) #the number of knots in the support
sp@der=list(matrix(rnorm(m*(k+1)),ncol=(k+1))) #the derivative matrix at random
sp1 = is.splinets(sp)[[2]] #the correction of 'der' matrices

support=list(matrix(c(5,12,17,29),ncol=2,byrow = TRUE))
sp=new("Splinets",knots=xi,smorder=k,supp=support) 
m=sum(sp@supp[[1]][,2]-sp@supp[[1]][,1]+1) #the number of knots in the support
sp@der=list(matrix(rnorm(m*(k+1)),ncol=(k+1))) #the derivative matrix at random
sp2 = is.splinets(sp)[[2]] 

spp = gather(sp1,sp2)

support=list(matrix(c(3,10,14,21,27,34),ncol=2,byrow = TRUE))
sp=new("Splinets",knots=xi,smorder=k,supp=support) 
m=sum(sp@supp[[1]][,2]-sp@supp[[1]][,1]+1) #the number of knots in the support
sp@der=list(matrix(rnorm(m*(k+1)),ncol=(k+1))) #the derivative matrix at random
sp3 = is.splinets(sp)[[2]] 

spp = gather(spp, sp3)

plot(spp)
gramian(spp) #the regular gramian matrix
spp2=subsample(spp,sample(1:3,size=3,rep=TRUE))
gramian(Sp=spp,Sp2=spp2) #cross-Gramian matrix


#-----------------------------------------#
#--------- Grammian for B-splines --------#
#-----------------------------------------#

n=25; xi=seq(0,1,by=1/(n+1)); k=2; 
Sp=splinet(xi) #B-splines and corresponding splinet

gramian(Sp$bs) #band grammian matrix for B-splines
gramian(Sp$os) #diagonal gramian matrix for the splinet
A=gramian(Sp=Sp$bs,Sp2=Sp$os) #cross-Gramian matrix, the coefficients of
                              #the decomposition of the B-splines

plot(Sp$bs)
plot(lincomb(Sp$os,A))

Splinets documentation built on March 7, 2023, 8:24 p.m.