Analysis of Longitudinal Data

Analysis of longitudinal data -an application to chronic angle
closure glaucoma
Pallavi Basu
Abhishek Pal Majumder
Anirban Basak
Priyam Biswas
May 29, 2007
Abstract
We examine longitudinal data of visual ?eld score and IOP from patients having chronic angle closure glaucoma. In
determinig a relationship between ?eld score and IOP , linear regression technique is used . Serious concerns can be raised
about the normality assumption. A Box-Cox transformation is hence applied.We try to analyze the assumption that each
sub?eld is equally a?ected by glaucoma .Resampling technique is used to estimate distribution of test statistic . Predicting
Progression was not feasible due to shortage of data.
1

---

Contents
1 Outlining situation and framing objectives.
4
1.1
Explaining the variables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4
1.2
Inclusion criterion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4
1.3
Categorization by glaucoma stage(by AGIS system) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4
1.4
Objectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4
2 Brief description of methods of analysis
4
3 Using this dataset
5
3.1
Dealing with missing data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5
3.2
One assumption that can’t be ignored here . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5
3.3
Handling of visual acuity
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5
4 Examining relationship between IOP and Visual ?eld score
5
4.1
Selection of response and Explanatory variables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5
4.2
Independence of left and right eye . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5
4.2.1
Preliminary Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5
4.2.2
Formulation of hypothesis and testing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5
4.2.3
Evaluation of cut-o? . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
6
4.2.4
Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
6
4.2.5
Interpretation of results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
6
4.3
Choice of model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7
4.4
Selection of structure of V0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7
4.4.1
Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7
4.4.2
The exponential correlation model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7
4.4.3
Justi?cation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7
4.5
Method of analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7
4.5.1
Restricted maximum likelihood estimation(REML) . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7
4.5.2
Box-Cox transformation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
8
4.5.3
Box-Cox transformation and REML . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
8
4.6
Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
8
4.6.1
Estimates of the parameters of the model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
8
4.6.2
Model adequacy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
8
4.7
Hypothesis testing
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
8
4.8
Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
10
4.9
Interpretation of results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
10
4.10 Testing between nonNONE categories and interpretation of results . . . . . . . . . . . . . . . . . . . . . . . .
10
4.11 An interesting observation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
10
5 To evaluate characteristic visual ?eld defect
11
5.1
De?ning baseline ?eld score . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
11
5.2
Methodology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
11
5.3
Category : MILD . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
11
5.3.1
Hypothesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
11
5.3.2
Testing procedure . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
11
5.3.3
Evaluation of cut-o? . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
11
5.3.4
Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
13
5.3.5
Interpretation of results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
13
5.4
Category : MODERATE
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
13
5.4.1
Hypothesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
13
5.4.2
Testing procedure . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
13
5.4.3
Evaluation of cut-o? . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
13
5.4.4
Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
13
5.4.5
Interpretation of results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
13
6 To evaluate Progression of visual ?eld damage
14
6.1
De?nition of Progression
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
14
6.2
Objectives and problem faced
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
14
6.3
Future scope . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
14
7 Dealing with missing data
15
7.1
Dropouts and intermittent missing values . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
15
7.2
Dealing with intermittent missing values . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
15
7.3
Methodology for dropouts
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
15
2

---

8 References
16
9 Acknowledgements
17
List of Figures
1
Scatter plot of left and right IOP at di?erent time points . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
6
2
Scatter plot of residual . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
9
3
Normal probability plot of residuals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
9
4
Empirical cdf for T1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
12
5
Empirical cdf for T2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
12
6
Empirical cdf for T . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
14
3

---

1
Outlining situation and framing objectives.
90 patients each having chronic angle closure glaucoma in one or both pairs of eyes were diagonised at 4 di?erent time points
within a time span of two years.The purpose of this project is to resolve out issues that are of help to medical experts.
1.1
Explaining the variables
1. Age of the patient at the ?rst time point of visit
2. Gender
3. Visual acuity
4. Intraocular pressure(IOP)
5. Field score(0-20)1
• Nasal (0-2)
• Superior hemi?eld(0-9)
• Inferior hemi?eld(0-9)
6. Kind of treatment provided upto that time point
• Drop
• Trab
• PI
• IOL
• Needling
1.2
Inclusion criterion
• New and follow up cases of chronic angle closure glaucoma with or without treatment
• Absence of any other major eye disease or any other kind of glaucoma
• Age group 30-70
1.3
Categorization by glaucoma stage(by AGIS system)
• None(score 0)
• Mild(score 1-5)
• Moderate(score 6-11)
• Severe(score 12-17)
• End-stage(score 18-20)
1.4
Objectives
1. To evaluate characteristic visual ?eld defect .
2. To assess the relationship between IOP and visual ?eld damage.
3. To evaluate Progression of visual ?eld damage.
2
Brief description of methods of analysis
• Under the null hypothesis it is expected that all the ?elds - Nasal , Superior hemi?eld and Inferior hemi?eld have same
e?ect. A suitable test statistic is used . The distribution of this statistic is estimated using large no. of resamples of
equal size from the original data2.A 100(1 ? 0.05)% con?dence interval is obtained. Depending on where the value
under null lies the null is rejected or not rejected.
• Linear regression model taking IOP as response and all other as explanatory variables is used.As the data of a particular
individual is correlated over time an autoregressive model of order 1 is selected. Treating the dataset as longitudinal
data, linear regression model is ?tted. Usual methods of analysis follows thereafter.
• Due to a very small no. of data, proper analysis of progression is not feasible. Many more follow ups of the same
datasets are required.
1 AGIS scoring system is universally used in this work
2 in literature called Bootstrap
4

---

3
Using this dataset
3.1
Dealing with missing data
Missing data will be categorised into dropouts and other.A separate section3 deals with this problem.However, not much
missing data - specially the dropout variety is present in this particular dataset.
3.2
One assumption that can’t be ignored here
Due to unavailability of data , equal time spaced interval for all subsequent visits and for all patients are assumed.Although,
this may seem a crude assumption taking into account that a medical expert asks a patient for a next visit at an interval
that a particular patient deserves.From this viewpoint, the assumption is partially justi?ed.
3.3
Handling of visual acuity
This has a special kind of data form.A chart4 is used to measure visual acuity.If a person has visual acuity 20/40, at 20 feet
from the chart that person can read letters that a person with 20/20 vision could read from 40 feet away.Since a linear model
is assumed , the visual acuity score is converted to fraction(20/40 ? 1/2).
4
Examining relationship between IOP and Visual ?eld score
4.1
Selection of response and Explanatory variables
IOP is treated as a response and all other variables as explanatory.
4.2
Independence of left and right eye
For each time point,a nonparametric approach to test independence is prefered.Kendall’s test for independence based on
signs is used 5.To avoid the di?culty due to ties the usual cut-o? is not used.Fixing an eye, the IOP values for the other eye
is permuted over di?erent patients6 to get an estimate of cut-o?.
4.2.1
Preliminary Analysis
The Correlation coe?cient between left eye and right eye at various time points are examined.Notice that the absolute values
decreases with time .
Time point 1:? = +0.3644
Time point 2:? = +0.3429
Time point 3:? = ?0.1037
Time point 4:? = ?0.0324
4.2.2
Formulation of hypothesis and testing
For each time point:
H0 ? ? = 0
and
H1 ? ? = 0
The Kendall sample correlation statistic for X and Y of n independent paired sample is
n?1
n
K =
Q((Xi, Yi), (Xj, Yj))
i=1 j=i+1
? 1 if (d?b)(c?a) > 0
where, Q((a,b),(c,d)) =? ?1 if (d ? b)(c ? a) < 0
? 0 if (d?b)(c?a) = 0
RejectH0 if K ? k?/2 or K ? k
where, k
is lower
?/2
?/2 is upper ?/2 tail probability of the null distribution of K and k?/2
?/2 tail probability of the null distribution of K.
3 refer section 8
4 technically called Snellen chart
5 for detailed theory refer Nonparametric statistical methods(Hollander and Wolfe)
6 termed ’Permutation distribution in literature
5

---

Time point 1
Time point 2
60
60
40
40
20
20
0
0
0
20
40
60
0
20
40
60
Time point 3
Time point 4
40
30
30
20
20
10
10
0
0
0
20
40
60
0
10
20
30
40
50
Figure 1: Scatter plot of left and right IOP at di?erent time points
4.2.3
Evaluation of cut-o?
A signi?cant proportion of the data has tied ranks.Due to this drawback usual cut-o? tables cannot be refered.Permutation
distribution was used to evaluate the cut-o?.Keeping IOP values of one of the eyes ?xed, the other IOP value was permuted
at random.The Kendall’s test statistic was recalculated using this dataset.This procedure is repeated for 10000 times. From
the emperical cdf of this statistic , 100(1 ? 0.05)% CI was constructed.The null is then rejected or not rejected accordingly
as the kendall’s statistic from the original dataset lies within or outside this CI.
4.2.4
Results
The results for the 4 time points are as follows :
• At time point 1 the null was rejected
• At time point 2 the null was rejected
• At time point 3 the null was not rejected
• At time point 4 the null was not rejected
4.2.5
Interpretation of results
From the sample correlaion coe?cient at various time points mentioned earlier it was already observed that its absolute
value decreases with time.Moreover, from the above result at the last two time points , the IOP values for the pair of eyes
of an inidividual is found to be uncorrelated. This hints to review the dataset.It is then observed that at the ?rst time point
almost 50% of the patients had no medical treatment before.This leads to consider the treatment nonuniformity among the
patients.This e?ect is further continued to the second time point.But,when at the third time point the nonuniformity of being
under medication or not disappears and hence for the third and the fourth time points the IOP values for the pair of eyes
are uncorrelated.Hence,if treatment is considered as an explanatory variable IOP values of pair of eyes of a patient can be
considered uncorrelated.
6

---

4.3
Choice of model
˜
YIOP = ?0 +
˜
Xage?1 +
˜
Xgender?2 +
˜
Xvisualacuity?3 +
˜
Xdrop?4 + ˜
Xiol?5 +
˜
Xtrab?6 + ˜
Xpi?7 +
˜
Xneedling?8 +
˜
XMILD?9 +
˜
XMODERAT E?10 +
˜
XSEV ERE?11 +
˜
XENDST AGE?12 + ˜
?
Since, the IOP is taken as a response variable and all others as explanatory variables, from the earlier conclusions the
left and right eyes of an individual are taken as independent experimental units.However, measurements of an unit over
di?erent time points cannot be taken uncorrelated.
GLM for longitudinal data treats y as a realization of a multivariate Gaussian random vector Y with
Y ? M V N (X?, ?2V )
where, V is a block diagonal matrix with nonzero 4 × 4 blocks V0 , each representing the variance matrix for the vector of
measurements on a single experimental unit.
4.4
Selection of structure of V0
4.4.1
Motivation
The sample time correlation matrix is:
1.00
0.51
0.40
0.32
0.51
1.00
0.50
0.44
0.40
0.50
1.00
0.51
0.32
0.44
0.51
1.00
Notice correlation between ?rst and second time point is almost same with that of second and third and also with third and
fourth time point.Moreover, correlation between ?rst and third time point is close to that of second and fourth time point.
4.4.2
The exponential correlation model
In this model V0 has jkth element, vjk = Cov(Yij, Yik) of the form vjk = ?2?abs(j?k)
Yij denotes the observation of ith experimental unit at jth time point.
A justi?cation of above model is to represent the random variables Yij as Yij = µij + Wij, i = 1, . . . , m ,j = 1, . . . , n,where
Wij = ?Wij?1 + Zij, and Zijs are mutually independent N(0,?2(1 ? ?2)) where m = no. of experimental units & n = no. of
time points for each unit .
4.4.3
Justi?cation
In the exponential correlation model the correlation between jth and kth time points of an individual depends on j and k only
through their absolute di?erence.As, the sample correlation matrix almost satis?es this property ,the exponential correlation
model is selected.
4.5
Method of analysis
4.5.1
Restricted maximum likelihood estimation(REML)
In the case of the GLM with dependent errors the REML estimtor is de?ned as a maximum likelihood estimator based on a
linearly transformed set of data Y ? = AY such that the distribution of Y ? does not depend on ? where
Y ? M V N (X?, ?2V )
. Calculation7shows that the REML estimator maximises the loglikelihood equation
L?(?2, V ) = ?0.5 log(det(?2V )) ? 0.5 log(det(??2X V ?1X)) ? 0.5(y ? X ˆ
?) ??2V ?1(y ? X ˆ
?)
Substituiting, ˆ
? and ˆ
?2 in the loglikelihood equation ,
L?(V0) = ?0.5m(n log(RSS(V0)) + log(det V0)) ? 0.5 log(det(X V ?1X))
where,
RSS(V0) = (y ? X ˆ
?(V0)) V ?1(y ? X ˆ
?(V0))
and
ˆ
?(V0) = (X V ?1X)?1X V ?1y
To solve V0 method of iteration is used.Subsequently, ˆ
? and ˆ
? are obtained.
7 Refer Analysis of longitudinal dataDiggle Chapter 4
7

---

4.5.2
Box-Cox transformation
Original data on IOP being integer valued it is wise to apply a Box-Cox transformation8 to ensure normality.
(y? ? 1)/? ?y??1
if ? = 0
y(?) =
?y ln y
if ? = 0
where ?y is the geometric mean of the response variable.Applying this transformation SSE(?) is calculated for di?erent values
of ? and that value of ? is chosen for which SSE(?) is minimum.
4.5.3
Box-Cox transformation and REML
There being no closed form solution of the REML log-likelihood equation and V0 being function of ? only, di?erent values of
? were used to evaluate the log-likelihood equation and ˆ
? is that which maximizes the log-likelihood equation.To implement
Box-Cox transformation in this set-up ?rst a ? is ?xed for which ˆ
? is evaluated and the corresponding SSE is obtained.Now,
varying over ? , ˆ
? is obtained for which SSE(?) is minimum.
4.6
Results
4.6.1
Estimates of the parameters of the model
ˆ
?=0.6
ˆ
?=0.36
29.22
?0.01
?1.42
?0.26
?3.06
5.21
ˆ˜
?= ?9.16
?1.84
?7.08
2.41
3.27
4.57
3.03
4.6.2
Model adequacy
• The scatter plot of the residuals(Figure 2 ) appears to be random which emphasizes that they do not exihibit any
de?nite pattern.
• The normal probability plot of the residuals (Figure 3 ) appears to be in a straight line indicating that the fact errors
are indeed normal,emphasizing normality assumption of the response is valid .
4.7
Hypothesis testing
H0 ? Q? = 0
,where Q is a full-rank q × p matrix for some q ? p. It can be deduced that
Q
ˆ
?REML ? M V N (Q?, Q
ˆ
RREMLQ )
,where
ˆ
RREML = ˆ
?2(X
ˆ
V ?1
X)?1
REM L
An appropriate test statistic for testing the hypothesis Q? = 0 would be
T =
ˆ
?REML Q (Q
ˆ
RREMLQ )?1Q
ˆ
?REML
and the approximate null sampling distribution of T is chi-squared on q degrees of freedom.
8 For details refer Design of experimentsMontogomery
8

---

Plot of residual
20
15
10
5
0
?5
?10
?15
?200
100
200
300
400
500
600
Figure 2: Scatter plot of residual
Normal Probability Plot
0.999
0.997
0.99
0.98
0.95
0.90
0.75
0.50
Probability
0.25
0.10
0.05
0.02
0.01
0.003
0.001
?15
?10
?5
0
5
10
15
Data
Figure 3: Normal probability plot of residuals
9

---

4.8
Results
It is required to ?nd out whether there is a ?eld category e?ect on the IOP.Each of the p-values listed below indicate the
result for :
H0 ?
ˆ
?correspondingcategory = 0
versus
H1 ?
ˆ
?correspondingcategory = 0
• P-value of
ˆ
?MILD = .029
• P-value of
ˆ
?MODERAT E = .006
• P-value of
ˆ
?SEV ERE ? 0
• P-value of
ˆ
?ENDST AGE = .02
4.9
Interpretation of results
A p-value less than 0.05 here, indicates that for any unit in that category the expected value of IOP is higher than that if
the unit was in the None category. Also, if somebody su?ers from glaucoma he/she is bound to have higher IOP than that
of a normal person.
4.10
Testing between nonNONE categories and interpretation of results
In similar lines, proper testing was done to check whether there is signi?cant di?erence among the four categories -
mild,moderate,severe and end stage.
E?ect of severe is found to be the most whereas e?ect of mild is found to be the least amongst the four categories.However,di?erence
of e?ect of moderate and that of end stage was not signi?cant.
Restating the two possible order of increasing e?ect of visual ?eld categories are,
mild < moderate < endstage < severe
or,
mild < endstage < moderate < severe
4.11
An interesting observation
Normal probability plot of the residuals showed that there were 8 outliers. Retrieving the original data it has been observed
that
• In cases where outliers have positive residuals trabeculectomy has been done just after the time point at which residual
is an outlier
• In cases where outliers have negative residuals trabeculectomy has been done just before the time point at which residual
is an outlier
This emphasizes the fact that trabeculectomy has an enormous e?ect in reducing the IOP of patients having glaucoma .
10

---