GLSME-package: Generalized Least Squares with Measurement Error

Description Details Author(s) References See Also Examples

Description

The code fits the general linear model with correlated data and observation error in both dependent and independent variables. The code fits the model

y = Dβ + r, r \sim N(0,V), V = σ^{2} T + V_{e} + Var[Uβ|D],

where y is a vector of observed response variables, D is an observed design matrix, β is a vector of regression parameters to be estimated, σ^{2}T is a matrix representing the true residual variance, V_{e} is a matrix of known measurement variance in the response variable, and Var[Uβ|D] is a matrix representing effects of measurement error in the predictor variables (see Hansen and Bartoszek 2012).

Details

Package: GLSME
Type: Package
Version: 1.0
Date: 2014-05-19
License: GPL (>= 2)
LazyLoad: yes

The code fits the general linear model with correlated data and observation error in both dependent and independent variables. The code fits the model

y = Dβ + r, r \sim N(0,V), V = σ^{2}T + V_{e} + Var[Uβ|D],

where y is a vector of observed response variables, D is an observed design matrix, β is a vector of regression parameters to be estimated, σ^{2}T is a matrix representing the true residual variance, V_{e} is a matrix of known measurement variance in the response variable, and Var[Uβ|D] is a matrix representing effects of measurement error in the predictor variables (see Hansen and Bartoszek 2012).

The estimation function is GLSME. It is an iterated (if the variance parameters are unknown) generalized least squares estimation procedure.

The motivation for the approach is that the observations and errors are correlated due to an underlying phylogeny but the program allows for any dependence structure.

In the mvSLOUCH package an alternative method of correcting for observation error is used. The error variance-covariance matrix enters the likelihood function by being added to the biological variance-covariance matrix.

Author(s)

Krzysztof Bartoszek Maintainer: <[email protected]>

References

Bartoszek, K. and Pienaar, J. and Mostad. P. and Andersson, S. and Hansen, T. F. (2012) A phylogenetic comparative method for studying multivariate adaptation. Journal of Theoretical Biology 314:204-215.

Hansen, T.F. (1997) Stabilizing selection and the comparative analysis of adaptation. Evolution 51:1341-1351.

Hansen, T.F. and Bartoszek, K. (2012) Interpreting the evolutionary regression: the interplay between observational and biological errors in phylogenetic comparative studies. Systematic Biology 61(3):413-425.

Hansen, T.F. and Pienaar, J. and Orzack, S.H. (2008) A comparative method for studying adaptation to randomly evolving environment. Evolution 62:1965-1977.

See Also

mvSLOUCH

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
library(mvSLOUCH)
library(ape)
library(ouch)
n<-5 ## number of species
apetree<-rtree(n)
phyltree<-ape2ouch(apetree) ##mvslouch requires ouch format
### Correct the names of the internal node labels.
phyltree@nodelabels[1:(phyltree@nnodes-phyltree@nterm)]<-
as.character(1:(phyltree@nnodes-phyltree@nterm))
### Define Brownian motion parameters to be able to simulate data under the Brownian motion model.
BMparameters<-list(vX0=matrix(0,nrow=2,ncol=1),Sxx=rbind(c(1,0),c(0.2,1)))
### Now simulate the data and remove the values corresponding to the internal nodes.
xydata<-simulBMProcPhylTree(phyltree,X0=BMparameters$vX0,Sigma=BMparameters$Sxx)
xydata<-xydata[(nrow(xydata)-n+1):nrow(xydata),]

x<-xydata[,1]
y<-xydata[,2]

yerror<-diag((rnorm(n,mean=0,sd=0.1))^2) #create error matrix
y<-rmvnorm(1,mean=y,sigma=yerror)[1,]
xerror<-diag((rnorm(n,mean=0,sd=0.1))^2) #create error matrix
x<-rmvnorm(1,mean=x,sigma=xerror)[1,]
GLSME(y=y, CenterPredictor=TRUE, D=cbind(rep(1, n), x), Vt=vcv(apetree), 
Ve=yerror, Vd=list("F",vcv(apetree)), Vu=list("F", xerror))

Example output

Loading required package: mvtnorm
Loading required package: corpcor
Loading required package: ouch
Loading required package: subplex
Loading required package: ape
Loading required package: numDeriv

Attaching package: 'mvSLOUCH'

The following object is masked from 'package:GLSME':

    .createCovariancematrix

[1] "Finished running iteration 1 of iterated GLS. Current estimate : "
                      [,1]
                 0.6606684
x                0.8748019
ResponseVariance 1.2211642
[1] "Finished running iteration 2 of iterated GLS. Current estimate : "
                 ResponseVariance
                        0.6609151
x                       0.8786942
ResponseVariance        1.1851081
[1] "Finished running iteration 3 of iterated GLS. Current estimate : "
                      [,1]
                 0.6609233
x                0.8782650
ResponseVariance 1.2365300
[1] "Finished running iteration 4 of iterated GLS. Current estimate : "
                 ResponseVariance
                        0.6609163
x                       0.8788820
ResponseVariance        1.2488974
$GLSestimate
       [,1]
  0.6609154
x 0.8790248

$errorGLSestim
       [,1]
  0.5783075
x 0.8129168

$BiasCorrectedGLSestimate
       [,1]
  0.6609249
x 0.8883596

$errorBiasCorrectedGLSestim
       [,1]
  0.5783075
x 0.8215496

$K
                x
  1 -0.0000106928
x 0  0.9894921226

$R2
          [,1]
[1,] 0.4983162

$R2BiasCorrectedModel
          [,1]
[1,] 0.4984191

$LogLikelihood
[1] -5.68026

$LogLikelihoodBiasCorrectedModel
[1] -5.680138

$PredictorVarianceConstantEstimate
                  x 
0.0000000 0.3831323 

$ResponseVarianceConstantEstimate
[1] 1.248897

GLSME documentation built on May 29, 2017, 9:31 a.m.