Revision 948ad5967298b4dc174f4d96511af60f38e9279e authored by Roger Koenker on 27 July 2019, 09:38:23 UTC, committed by cran-robot on 27 July 2019, 09:38:23 UTC
1 parent 2b3a85a
rls.r
#This is a simple recursive least squares routine Reference: Harvey TSM p. 100
#
subroutine rls(n,p,x,y,b,A,Ax)
integer n,p
double precision x(p,n),y(n),b(p,n),A(p,p),Ax(p)
double precision zero,one,mone,f,r,ddot
data one/1.d0/
data mone/-1.d0/
data zero/0.d0/
#
#On input:
#
# A = crossprod(x[1:p,1:p))^{-1}
# b(,p) = A crossprod(x[1:p,1:p],y[1:p])
#
do i = (p+1),n {
call dgemv('N',p,p,one,A,p,x(1,i),1,zero,Ax,1)
f = one + ddot(p,x(1,i),1,Ax,1)
r = (y(i)-ddot(p,x(1,i),1,b(1,i-1),1))/f
call daxpy(p,one,b(1,i-1),1,b(1,i),1)
call daxpy(p,r,Ax,1,b(1,i),1)
call dger(p,p,mone/f,Ax,1,Ax,1,A,p)
}
return
end
Computing file changes ...