# R/compute.CTRinvC.r In SemiPar: Semiparametic Regression

```########## compute.CTRinvC.S #######################################

# Computes the inverse residual covariance matrix R
# for linear mixed models

# inputs: individual = id's for the "individiduals" on whom
#         the repeated mesaures are taken. This vector has
#         the same length as the response and is assumed to
#         be sorted into homogenous groupings. i.e.:
#         1 1 1 2 2 3 3 3 3 3 is OK
#         1 2 3 1 2 3 1 3 3 3 is not
#         rho = AR(1) parameter
#
#
# output: t(C) R.unscaled.inv C
#         C = [X Z]
#
# Last changed: 10/20/00

compute.CTRinvC <- function(X,Z,individual,rho=0)
{
C <- cbind(X,Z)
CTRinvC <- 0*t(C)%*%C

# Calculate indices (matrix location) and time indices for each individiual

if (length(individual)>1)
{
end.indiv <- (1:length(individual))[diff(individual)!=0]
strt.indiv <- c(1,end.indiv+1)
end.indiv <- c(end.indiv,length(individual))
n.indiv <- length(end.indiv)
indiv.inds <- NULL
indiv.times <- NULL

for (i in 1:n.indiv)
{
indiv.inds <- (strt.indiv[i]:end.indiv[i])
indiv.times <- (1:length(indiv.inds))
R.inv <-
solve(rho^abs(outer(indiv.times,indiv.times,"-")))
if (ncol(R.inv)>1)
{
CTRinvC <- CTRinvC +
t(C[indiv.inds,]) %*% R.inv %*% C[indiv.inds,]
}
else
{
CTRinvC <- CTRinvC +
as.numeric(R.inv)*outer(C[indiv.inds,],C[indiv.inds,])

}
}
}
else
{
CTRinvC <- t(C)%*%C
}

return(CTRinvC)
}

########## end compute.CTRinvC.S #######################################
```

## Try the SemiPar package in your browser

Any scripts or data that you put into this service are public.

SemiPar documentation built on May 2, 2019, 5:42 a.m.