Hello, dear friend, you can consult us at any time if you have any questions, add WeChat: daixieit

MATH5885 Longitudinal Data Analysis

Week 2

Exercise 1

Read in the data, creating a new variable for log(FEV1/Height), and have a look at the first 5 rows:

>  fev  <- read .table(" . ./ . ./data/fev1 .txt", header = TRUE)

>  fev$lfev1h  <- fev$lfev1 - log(fev$height)

>  head(fev,  n  =  5)

subjectid height         age  iheight      iage      lfev1         lfev1h

1

1

1.20    9.3415

1.2  9.3415  0.21511  0.03278844

2

1

1.28  10.3929

1.2  9.3415  0.37156  0.12469992

3

1

1.33  11.4524

1.2  9.3415  0.48858  0.20340106

4

1

1.42  12.4600

1.2  9.3415  0.75142  0.40076313

5

1

1.48  13.4182

1.2  9.3415  0.83291  0.44086791

Question 1

Extract the residuals from the cubic spline fit:

>  fev$resid  <- resid(smooth .spline(fev$age, fev$lfev1h))

We can produce a spaghetti plot of the residuals:

> plot(fev$age[fev$subjectid  ==  1],  fev$resid[fev$subjectid  ==  1],

+

+

type="l",  lty=1,  col=1,  xlim  = range(fev$age),  ylim  = range(fev$resid), xlab="Age  (years)",  ylab="log(FEV1/Height)  residual")

> points(fev$age[fev$subjectid  ==  1],  fev$resid[fev$subjectid  ==  1], +               type="p", pch=1,  col  =  1,  cex  =  0 .6)

>  for  (i  in  2:max(fev$subjectid))  {

+      lines(fev$age[fev$subjectid  ==  i],  fev$resid[fev$subjectid  ==  i],  +                 type  =  "l",  lty  =  1,  col  =  i)                                                             +     points(fev$age[fev$subjectid  ==  i],  fev$resid[fev$subjectid  ==  i], +                   type  =  "p", pch  =  1,  col  =  i,  cex  =  0 .6)

+  }

>  text(fev$age[fev$subjectid  ==  197]  +  0 .8,

+           fev$resid[fev$subjectid  ==  197],  "id  197")

>  title("Spaghetti  Plot:  Residuals  log(FEV1/Height) minus  cubic  spline")

Spaghetti Plot: Residuals log(FEV1/Height) minus cubic splin

id 197

18

Age (years)

Question 2

First round the ages:

>  fev$ageround  <- round(fev$age)

In order to use the pairs function, we need to get the data into wide’ format, with one row per individual and a column for each age. We can do this using the reshape function:

> panel .cor  <- function(x, y, digits=2, prefix="", cex .cor) + {

+      usr  <- par("usr"); on .exit(par(usr))

+     par(usr  =  c(0,  1,  0,  1))

+     r  <- abs(cor(x, y,use="pairwise .complete .obs")) + txt <- format(c(r, 0 . 123456789), digits=digits)[1]

+      txt  <- paste(prefix, txt, sep="")

+      if(missing(cex .cor))  cex  <- 0 .8/strwidth(txt)

+      text(0 .5,  0 .5,  txt,  cex  =  cex  *  r)

+  }

>  fev .resid .wide .reshape  <- reshape (

+      fev[,c("subjectid",  "resid",  "ageround")],

+      direction  =  "wide",

+      v .names  =  "resid",

+      timevar  =  "ageround",

+      idvar  =  "subjectid")

>  head(fev .resid .wide .reshape)

subjectid         resid .9       resid .10       resid .11       resid .12       resid .13

1

8

16

1  -0 .18275144  -0 .15477209  -0 .15288553  -0 .03693050  -0 .07181568

2                   NA                   NA                   NA                   NA  -0 .10286747 3    0 .14248385    0 .19768037    0 .25395274                   NA  -0 .04114325


25                 4  -0 .13857584  -0 .06967834  -0 .06919457  -0 .10017548  -0 .18301857 35                 5                   NA                   NA                   NA                   NA                   NA 42                 6    0.06273254    0.11656534    0.06271131    0.04941064    0.00250609

resid .15       resid .16           resid .7         resid .8         resid .14         resid .17

1    -0 .13769410  -0 .18744768                     NA                   NA                     NA                     NA 8     0 .01590316    0 .07632892    0 .167993829    0 .06849432  -0 .096306567  -0 .168731997 16    0 .02461597  -0 .02089331    0 .219756542    0 .04612666  -0 .004626005                     NA 25  -0 .03473378                   NA  -0 .114311932  -0 .12945013    0 .008968504                     NA 35  -0 .16260641  -0 .11006280  -0 .084387049  -0 .12803007  -0 .057921453  -0 .155794884 42    0.02490978    0.06420085    0.001955933    0.08270616    0.125522304    0.004705282

resid .18  resid .6  resid .19

1                     NA           NA             NA 8     0 .03833939           NA             NA

16                   NA           NA             NA 25    0 .05364655           NA             NA 35  -0 .25670811           NA             NA

42                   NA           NA             NA

> pairs(fev .resid .wide .reshape[,- 1],  upper .panel  = panel .cor  )


−0.3                  −0.2                −0.4                −0.3                −0.4                −0.3                −0.20