Dokument Nr. K27453 aus den Wissensarchiven des GRIN Verlags

Kommentare und Fragen richten Sie bitte an:

E-Mail: info@grin.com
http://www.grin.com

von

Tony Wragg

Unsere Texte finden Sie in unseren Online-Archiven

www.GRIN.com

www.hausarbeiten.de

und auf angeschlossenen Partnersites.

---

Univariate & Multivariate Methods for the

Analysis of Repeated Measures Data.

Anthony J. Wragg

A thesis submitted in partial fulfilment of the requirements for the degree of Master of

Applied Science (Statistics and Operations Research).

Department of Statistics and Operations Research

Royal Melbourne Institute of Technology

December 1999

---

2

Declaration

The work contained in this thesis has not been submitted previously, in whole or in

part, in respect of any academic award.

To the best of my knowledge and belief, this thesis contains no material previously

published or written by any other person except whee due reference is given.

Anthony J. Wragg

31st December 1999

---

3

Acknowledgments

I would like to thank my thesis supervisor Ms. Kaye E. Marion and Associate

Professor Panlop Zeephongsekul for their encouragement and guidance in the writing

of this minor thesis and throughout the rest of the course.

---

4

Abstract

This thesis considers both univariate and multivariate approaches to the analysis of a

set of repeated-measures data. Since repeated measures on the same subject are

correlated over time, the usual analysis of variance assumption of independence is

often violated. The models in this thesis demonstrate different approaches to the

analysis of repeated-measures data, and highlight their advantages and disadvantages.

Milk from two groups of lactating cows, one group vaccinated, the other not, was

analysed every month after calving for eight months in order to measure the amount

of bacteria in the milk. The primary goal of the experiment was to determine if a

vaccine developed by the Royal Melbourne Institute of Technology′s Biology

Department led to a significant decrease in mean bacteria production per litre of milk

produced compared to the control group.

A univariate model suitable for repeated measures data was initially tried, with mean

bacteria production in the treatment group not significantly different from the control

group (

p

< 0.68).

The multivariate approach to repeated measures, profile analysis, yielded similar

results for treatment effects (

p

< 0.68), while meeting the necessary assumptions for

multivariate analysis.

Finally, a generalised multivariate analysis of variance was carried out in order to fit

polynomial growth curves for both the control and the vaccinated groups and to test if

the growth curves were equal for the two groups. It was found that a slope-intercept

model was adequate to describe both growth curves and that the growth curve for the

treatment group did not differ significantly from that of the control group (

p

< 0.11).

---

5

Table of Contents

1. INTRODUCTION 6

2. EXPLORING THE DATA 8

3. TIME BY TIME ANOVA 11

4. UNIVARIATE APPROACH TO REPEATED MEASURES 13

4.1 Repeated Measures ANOVA 13

4.2 Testing For Compound Symmetry 21

5. MULTIVARIATE APPROACH TO REPEATED MEASURES 26

5.1 Profile Analysis 28

5.2 Assumptions of Profile Analysis 29

5.3 Testing For Multivariate Normality 30

5.4 Testing For The Equality of Covariance Matrices 33

5.5 Hypothesis Tests for Profile Analysis 35

6. THE GENERALISED MULTIVARIATE ANALYSIS OF VARIANCE 42

6.1 Growth Curves 42

6.2 Hypothesis Tests For Growth Curves 48

6.3 Testing Polynomial Adequacy 49

6.4 Testing For The Equality of Growth Curves 53

7. CONCLUSION 59

REFERENCES 61

APPENDIX A 63

---

6

1. Introduction

Milk from two groups of lactating cows, one group vaccinated, the other not, was

analysed every month after calving for eight months in order to measure the amount

of bacteria in the milk. The primary goal of the experiment was to determine if a

vaccine developed by the Royal Melbourne Institute of Technology′s Biology

Department led to a significant decrease in mean bacteria production compared to the

control group.

Experiments such as this fit into the family of designs known in the literature as

repeated measures data

,

longitudinal models

, or

growth curves

. Data from these

models generally arise whenever more than two observations of the same variable are

made on an individual subject or experimental unit. These models are especially

common in biology, agriculture, and medicine and most often occur when

observations on a group of subjects are repeated over a period of time.

Repeated measures data such as this require somewhat different statistical treatment

than normal because the observations are not independent. This lack of independence

lies at the core of repeated measures analysis, and is what differentiates it from the

more commonly used statistical analyses. The implications of a lack of independence

within a subject′s responses are serious, as will be explained in Chapter 4. For the

data collected by the RMIT biology department, this implies that the amount of

bacteria found in a litre of a cow′s milk will, at any given time, be correlated with the

amount of bacteria found in that cow′s milk at subsequent or preceding times. In

addition, the correlation between the amount of bacteria produced at different times

also tends to be stronger the shorter the time interval. In other words, the amount of

bacteria produced per litre of milk is more dependent on the amount of bacteria

produced one month ago than the amount of bacteria produced five months ago.

Correlation between observations is usually present in these types of experiments.

Nearby plots in a field trial are usually more similar than plots further apart.

---

7

When applying different levels of a factor, the effects of this correlation are generally

overcome by randomisation ­ the levels are randomly allocated to the experimental

subjects. Randomisation ensures that in the long run there is no correlation between

the factor levels, so that observations with any given factor level are not more similar

to some factor levels than to others. Since time is treated as a factor with the eight

months considered the eight levels of time, randomisation is impossible ­ the

observations must follow their natural sequence. Thus, it is not possible to randomise

the order of monthly observations: they must follow the sequence month 2, month 3,

month 4... As a result of the lack of randomisation, the means of two milk samples

taken a month apart, for example, tend to be more highly correlated than those taken 6

months apart. As a consequence, the precision of the difference tends to drop as the

time interval increases, nullifying the use of a single standard error of difference for

the time factor.

This variation in correlation between levels of the time factor means that it is

inappropriate to analyse the data as if time was a randomised factor. This thesis

covers some of the most commonly used techniques. Frequently, there is no single

best approach to analysis. It depends on what questions need to be answered. Often, it

is useful to use two or more approaches with the same data.

The analysis carried out on a set of repeated measures data is determined largely by

the questions the researcher wants answered. For this study, the most important

question is that concerning the vaccine: is there a significant difference between the

mean number of cells found in a litre of milk produced by cows in the treatment group

compared to the control group, over time? Secondary questions may include such

things as the change in cell production over the months irrespective of group

membership. The first question when studying repeated measures, or in fact any, data,

should not be how to analyse the data, but what is the experimenter is interested in

finding out (Lindsey, 1993). Once this is known, together with knowledge of the

techniques available, the selection of an appropriate technique becomes much easier.

---

8

2. Exploring The Data

The biology department of RMIT has developed a vaccine which is thought to reduce

the number of cells of mastitis, hereto known as `cells′ found in a cow′s milk. The

vaccine was tested on a randomly chosen sample of 23 cows, while a randomly

chosen sample of 18 cows was used as the control. Readings of the cell count of each

cow′s milk were taken at 2, 3, 4, 5, 6, 7, 8, and 9 months after calving.

One of the attractive features of repeated measures data is that they can be displayed

in a graphical plot which is readily interpretable, without requiring a great effort and

little training is required to interpret the plots (Lindsey, 1993). Data plotting is

essential in order to get some feeling for what patterns are present in the data, whether

expected trends have occurred, what unexpected features are apparent and what

questions deserve analytical consideration. It should always precede detailed analysis

of the data. Figure 1 (a)-(d) shows various plots of the two groups.

Figure 1.

CONTROL

TREATMENT

3500

3500

3000

3000

2500

2500

2000

2000

1500

1500

1000

1000

500

500

0

0

2

3

4

5

6

7

8

9

2

3

4

5

6

7

8

9

MONTH

MONTH

(a)

(b)

3500

3500

3000

3000

2500

2500

2000

2000

LLS

LLS

1500

E

1500

E

C

C

1000

1000

500

500

0

0

2

3

4

5

6

7

8

9

2

3

4

5

6

7

8

9

MONTH

MONTH

(c)

(d)

---

9

The skew towards higher values of cells can be seen over most of the time periods,

although some months are worse than others. The boxplots, where the boxes contain

50% of the data, tend to vary over time in the control group. The pattern of variation

is roughly similar for the treatment group. The treatment group appears to have

comparatively smaller variances over time than the control group. If it were not for

the larger variance at month 3, the variance of the control group would be increasing

with time. The variance of the treatment group looks as if it decreases slightly over

time. Plots of the individual cow′s results show that the control group cows cell

production increases linearly over time, although the trend is not patent. The treatment

group does not show any apparent increase over time. There appears to be little

evidence of a quadratic or cubic growth curve from the plots.

The outliers were noted and checked for accuracy with the RMIT Biology department

to make sure there were no transcription errors or similar problems. The outliers were

all legitimate observations.

The sample means of each group are plotted in Figure 2(a). Plotting the medians,

Figure 2(b) as well as the means allows one to look at the data from a slightly

different perspective, one that is resistant to outliers. Since the observations for any

given month are generally skewed, the medians are a useful adjunct. The control

group′s mean response is not as stable as the treatment group′s.

Figure 2.

MEANS

MEANS

1500

CELLS

1500

CELLS

1000

1000

500

500

0

0

2

3

4

5

6

7

8

9

2

3

4

5

6

7

8

9

MONTH

MONTH

(a)

(a)

Control - - - - - Treatment

---

10

Fitting ordinary linear regression equations to the data gives an indication of the linear

trend for each group. Figure 3 (a) shows the control group has a positive slope which

highlights the increase in cells over time, while Figure 3(b) shows the treatment

group′s slope is negative and flatter. At first it might appear that a possible model for

these observations is the general linear model. The problem is, like ANOVA, that the

assumptions of linear regression require independence of the variables. In chapter 6

growth curves will be fitted to the data using a multivariate approach which does not

have the restrictive assumption of independence.

Figure 3.

Regression Plot

Regression Plot

Y = 280.095 + 48.0218X

Y = 528.168 - 2.35559X

R-Sq = 3.9 %

R-Sq = 0.0 %

3500

2500

3000

2000

2500

S

1500

2000

S

L

L

L

L

1500

CE

1000

CE

1000

500

500

0

0

2

3

4

5

6

7

8

9

2

3

4

5

6

7

8

9

MONTH

MONTH

(a) (b)

Although plotting the data is imperative with any analysis, it does not allow the

analyst to make anything more than general statements summarizing the apparent

behaviour of the subjects being studied. What is really needed is to be able to quantify

the responses and formally test the research questions given in chapter 1 as

hypotheses. Chapter 3 describes the most simple method for doing this.

---

11

3. Time by Time ANOVA

One of the simplest forms of longitudinal analysis is a time-by-time ANOVA. It

consists of

p

separate analyses, on for each subset of data corresponding to each time

of observation

t

. For more than two groups each analysis is a conventional ANOVA,

however since there are only two groups being compared in this study, the ANOVA

reduces to a two-sample

t

-test of

H

0

:

µ

control

= µ

treatment

at each of the

p =

8 times of

measurement. Table 1 shows the time-by-time ANOVA results for the data

Table 1.

Month 2

3

4

5

6

7

8

9

t

-0.80 0.14 -1.12 -1.56 -0.45 -0.47 1.02 1.21

p-value 0.43 0.89

0.27 0.13 0.65 0.64 0.32

0.24

The time-by-time analysis indicates that mean cells count does not differ between the

control and treatment groups in any of the 8 months. This suggests that the two mean

response profiles are alike. Month 5 has a large

t

test statistic (

t

= -1.56), but not

enough to be significant.

A time-by-time ANOVA is reasonably clear and uncomplicated, however Diggle,

Liang and Zeger (1994) point out its two major limitations. Firstly, it cannot address

questions concerning treatment effects which relate to the longitudinal growth of the

mean response curves, i.e. the growth rates between successive months. Secondly, the

inferences made within each of the

p

separate tests are not independent of each other,

nor is it clear how they should be combined. For example, a succession of marginally

significant group mean differences may be compelling with weakly correlated data,

but much less so if there are strong correlations between successive observations on

each cow.

---

12

The principal virtue of the time-by-time ANOVA approach to longitudinal studies is

its simplicity. The computational operations are elementary and the approach uses

familiar procedures in finding a solution to the problem.

In summary, whilst the time-by-time ANOVA may be useful in particular

circumstances, Diggle, Liang and Zeger (1994), do not recommend it as a viable

approach to longitudinal data analysis.

---

13

4. Univariate Approach to Repeated Measures

A more sophisticated approach than the time-by-time

t

-tests is a

repeated measures

analysis of variance (ANOVA). In contrast to the time-by-time approach, Diggle,

Liang and Zeger (1994), regard it as a first attempt to provide a single analysis of a

complete longitudinal data set.

4.1 Repeated Measures ANOVA

Experiments utilising repeated measures designs differ from the usual ANOVA

models in that the levels of time cannot be randomly assigned to one or more the

experimental units in the experiment. In this case the levels of time cannot be

assigned at random to the time intervals, and thus the usual ANOVA models may not

be valid. Because of the non-random assignment of time, the errors corresponding to

the respective experimental units may have a covariance matrix which does not

conform to those for which the usual ANOVA analysis are valid.

The inherent dependence that is associated with repeated measures data introduces

extra complications into the analysis. Unfortunately, the simplifying properties arising

from data which are independently and identically distributed can no longer be relied

upon. To yield conclusions which are valid the analyst must take into account the

possible dependence within subjects. Fortunately, Diggle, Liang and Zeger (1994) and

Vonesh and Chinchilli (1997) have outlined methods which modify the problem so

that independence based methods like ANOVA can be used.

Superficially, the

N

cows by

p

months structure of the data resembles that of a

randomised block or split plot design, so there is a temptation to carry out a standard

two factor group × month ANOVA. Using the standard ANOVA approach to this

problem presents problems for the unwary. Employing a standard ANOVA model

would regard the control and treatment groups as a factor on two levels, and

more

importantly

, it would regard time as a factor on

p

levels. One of the difficulties with

this approach is that the allocation of times to the

p

observations within each cow

cannot be randomised.

---

14

As was mentioned in the introduction, randomisation of the various levels of a factor

is an essential requirement of ANOVA. Usually treatment factors are randomised

within blocks. This is assumed in ANOVA. For example, the various treatments in

Block 1 can be randomised.

BLOCK 1

D

F

A

E

B

C

H

G

COW 1

Month 4 Month 2 Month 9 Month 3 Month 5 Month 7 Month 6 Month 8

The times, if they are regarded as factor levels, cannot be randomised as above, an

essential requirement of ANOVA. They must follow their natural sequence. The first

measurement is the first measurement, it cannot be taken third. In general, these extra

complications mean that the simple univariate ANOVA

F

-tests will no longer be

valid.

An approach that takes advantage of the fact that the

p

= 8 measurements on each cow

are repeated observations on the same response variable, namely, cells, is given in

Vonesh and Chinchilli (1997). The additive model is

y

= µ + +

+ +

( ) +

ijk

j

i

(

j

)

k

jk

ijk

---

15

where

i

= 1,...,

nj

is an index for cow

i

within group

j

,

j

= 1,...,

r

is an index for group membership,

k

= 1,...,

p

is an index for levels of time,

yijk

is the response measurement at time

k

for cow

i

within group

j

,

µ is the overall mean,

j

is the added effect for group

j

,

k

is the added effect for time

k

,

i

(

j

) is the random effect due to cow

i

within group

j

,

()

jk

is the added effect for the group

j

× time

k

interaction, and

ijk

is the random effect on time

k

for cow

i

within group

j

.

In the general model for this experiment, 18 cows were randomly assigned to the

control group and 23 cows were randomly assigned to the treatment group, and each

cow′s milk was tested on eight occasions (months).

Geisser (1980) gives the form of the ANOVA table

---

16

Anova Table

_____________________________________________________________________

Source df S.S. F

_____________________________________________________________________

SSM

Month

p

- 1

SSM

(

N

-

r

)

SSE

(

N

-

r

)

SSG

Group

r

­ 1

SSG

(

r

- )

1 SSC

(

G

)

Cows(within Group) (

N

-

r

)

SSC

(

G

)

(

N

-

r

)

SSG

×

SSM

Group × Month (

p

- 1)(

r

- 1)

SSG

×

SSM

(

r

-

1

)

SSE

Error (

p

- 1)(

N

-

r

)

SSE

_____________________________________________________________________

Total

Np

­ 1

SST

_____________________________________________________________________

Note that in this model there are two error terms, where

represents the random

i

(

j

)

error due to cow

i

within group

j

. In addition, there is the assumption that the

′s

i

(

j

)

and the ′s are independent with

ijk

.

iid

.

N

( ,

0

2

) and .

i

.

id

.

N

( ,

0

2

)

ijk

i

(

j

)

For this approach to be appropriate, Milliken and Johnson (1984) claim that the

variance structure of = Cov( ) must satisfy the assumption of compound

ij

symmetry.

---

17

A covariance matrix

is of compound symmetry form if it can be expressed as

1

L

1

L

= 2

1

L

M

M M

M

M

1

The compound symmetry condition implies that the random variables are equally

correlated and have equal variances. In other words, the variances of the differences

between pairs of errors, such as - are equal for all

k

and

k

,

k

k

. The

ijk

ijk

variance structure is also called the uniform-variance, equi-variance or the equi-

correlation structure.

Assuming that the univariate repeated measures approach is appropriate, i.e. has the

compound symmetry structure described above, the

F

tests for the hypothesis of

parallelism (no group × month interaction), coincidence (no differences between

groups), and constancy (no differences among months) are more powerful than the

corresponding multivariate tests where no structure is assumed for (Vonesh &

Chinchilli, 1997).

Tests exist for compound symmetry. Firstly, the above model will be applied to the

raw data and an analysis of the residuals for the usual assumptions done before testing

the variance structure of the residuals for compound symmetry.

The results from carrying out the analysis on Minitab are

---

18

General Linear Model

Factor Type Levels Values

month fixed 8 2 3 4 5 6 7 8 9

group fixed 2 0 1

cow(group) random 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14

15 16 17 18 19 20 21 22 23 24 25 26 27

28 29 30 31 32 33 34 35 36 37 38 39 40 41

Analysis of Variance for cells, using Adjusted SS for Tests

Source DF Seq SS Adj SS Adj MS F P

group 1 67952 67952 67952 0.10 0.758

cow(group) 39 27522342 27522342 705701 4.94 0.000

month 7 2383923 2874711 410673 2.88 0.007

month*group 7 2178344 2178344 311192 2.18 0.036

Error 273 38992902 38992902 142831

Total 327 71145463

The ANOVA table shows a significant effect for the group × month interaction, as

well as month alone. Before any conclusions can be made, an analysis of the residuals

needs to be made. Figure 4 (a) shows the residuals plotted against the fitted values

from the above model. There is evidence of increasing variance as the predicted

values increase. There also appear to be many outliers with larger than expected

positive values. Also alarming is the obvious lack of normality of the residuals in

Figure 4(b).

---

19

Figure 4.

Residuals Versus the Fitted Values

Normal Probability Plot of the Residuals

3000

3

2

2000

e

1

c

or

S

1000

s

i

dual

0

al

e

R

m

or

N

-1

0

-2

-1000

-3

0

1000

2000

-1000

0

1000

2000

3000

Fitted Value

Residual

(a)

(b)

Since the residuals from the raw data are highly skewed as well as exhibiting

increased variance with higher values of fits, and a strong departure from normality, a

transformation was applied. The raw data were transformed by applying log

e

to the

cell counts. The transformed data are shown in Figure 5 (a)-(d).

Figure 5.

CONTROL

TREATMENT

9

9

8

8

)

7

S

)

7

L

L

6

LLS

E

6

E

(

C

ln

5

5

l

n

(

C

4

4

3

3

2

2

2

3

4

5

6

7

8

9

2

3

4

5

6

7

8

9

MONTH

MONTH

(a)

(b)

9

9

8

8

7

7

)

)

S

S

6

6

L

L

L

L

E

E

5

5

(

C

(

C

ln

ln 4

4

3

3

2

2

2

3

4

5

6

7

8

9

2

3

4

5

6

7

8

9

MONTH

MONTH

(c)

---

20

The log

e

transformation has had the effect of making the data more symmetrical, in

that it does not have the highly skewed values previously seen with the raw data. Most

of the outliers have also disappeared. Note that the treatment group still has a smaller

variance over time than the control group.

The univariate ANOVA was then re-run with the transformed data. The transformed

data produced the following results for the repeated measures ANOVA .

General Linear Model

Analysis of Variance for ln(cells), using Adjusted SS for Tests

Source DF Seq SS Adj SS Adj MS F P

group 1 0.3232 0.3232 0.3232 0.17 0.684

cow(group) 39 74.9833 74.9833 1.9226 5.44 0.000

month 7 15.4230 16.6787 2.3827 6.75 0.000

month*group 7 3.3003 3.3003 0.4715 1.33 0.234

Error 273 96.4203 96.4203 0.3532

Total 327 190.4500

The transformed data `ln(cells)′ produced a different set of test statistics than the raw

data. Most importantly, the month × group interaction is no longer significant (p =

0.234). This indicates that the log

e

of the number of cells found in a cow′s milk over

time does not depend on whether the cow was vaccinated or not. The main effect for

month is still highly significant (p = 0.000), indicating that at least two of the eight

months have unequal means. As with the raw data, there was no significant

difference in the mean number of cells produced between the two groups, ignoring

time, (p=0.684). The residual plots of the transformed data are seen in Figure 6(a) ­

(d).

---

21

Figure 6.

Residuals Versus the Fitted Values

Normal Probability Plot of the Residuals

3

2

2

1

e

1

c

or

0

S

s

i

dual

0

al

e

R

m

-1

or

N

-1

-2

-2

-3

-3

4.5

5.5

6.5

7.5

-3

-2

-1

0

1

2

Fitted Value

Residual

(a)

(b)

Residuals Versus month

Residuals Versus group

2

2

1

1

0

0

s

i

dual

s

i

dual

e

e

R

R

-1

-1

-2

-2

-3

-3

2

3

4

5

6

7

8

9

CONTROL

TREATMENT

month

(d)

(c)

The residuals for the transformed data do not exhibit the skewness of the raw data, nor

do the residuals increase in variance as the fitted values increase. The normal

probability plot also shows marked improvement, although the tails are less than

satisfactory. Many transformations were tried, but with little impact on the tails. The

variance of the residuals over both month and group appear reasonably uniform.

4.2 Testing For Compound Symmetry

Some software packages, SPSS, for example, test for compound symmetry, while

others, like Minitab, do not. Since Minitab was used for this analysis, testing for

compound symmetry was carried out using the method outlined by Milliken and

Johnson (1984). The compound symmetry test can be applied to test the hypothesis

that is compound symmetric, i.e. the variance-covariance matrix of the residuals is

approximately

---

22

1

L

1

L

= 2

1

L

M

M M

M

M

1

This is a likelihood ratio test and it rejects

H

0:

is compound symmetric, for

1 where the range of is

1 1

p

-1

An estimate of can be derived from the residuals by using

2

2

^

= [

p

(

s

)

ii

(

p

- )

1

p

p

s

2

i

=1

j

=

ij

1

]

where

sij

is the

ij

th element of the matrix

S

and

s

is the mean of the diagonal

ii

elements of

S

. The residuals for cow

i

in group

j

are represented as

^

ij

1

^

ij

2

^ = ^

ij

ij

3

M

^

ijp

The variance-covariance matrix of the residuals can be calculated as

r

n j

^ ^

j

=

1

i

1

ijk ijk

′

s

=

=

kk

′

N

-1

---

23

The residual variance-covariance matrix for the transformed data is

656

.

0

072

.

0

- 054

.

0

- 073

.

0

- 180

.

0

- 161

.

0

- 131

.

0

- 126

.

0

072

.

0

403

.

0

064

.

0

- 080

.

0

- 145

.

0

- 072

.

0

- 123

.

0

- 119

.

0

- 054

.

0

064

.

0

258

.

0

011

.

0

- 085

.

0

- 016

.

0

- 084

.

0

- 092

.

0

- 073

.

0

- 080

.

0

011

.

0

257

.

0

013

.

0

- 043

.

0

- 046

.

0

-

038

.

0

S

=

- 180

.

0

- 145

.

0

- 085

.

0

013

.

0

241

.

0

017

.

0

109

.

0

029

.

0

- 161

.

0

- 072

.

0

- 016

.

0

- 043

.

0

017

.

0

175

.

0

028

.

0

073

.

0

- 131

.

0

- 123

.

0

- 084

.

0

-

046

.

0

109

.

0

028

.

0

195

.

0

051

.

0

- 126

.

0

- 119

.

0

- 092

.

0

- 038

.

0

029

.

0

073

.

0

051

.

0

222

.

0

The test statistic for the compound symmetry test is

Q =

(1-

A

)

M

Where

Q

is approximately distributed as a 2

random variable with

f

degrees of

freedom where

M

= -

v

log

e

(^ )

p

(

p

+ )

1 2(2

p

- )

3

A

=

6

v

(

p

- )(

1

2

p

+

p

- 4)

2

p

+

p

- 4

f

=

2

and

v

is the degrees of freedom associated with the elements in . To ensure division

by zero does not occur in calculating the formula for

A

, the condition

p

2 must

hold. If

Q

> 2

then the null hypothesis is rejected, and the resulting

F

tests will

,

f

no longer be valid. In the case where

Q

is significant, then the usual procedure is to

multiply the degrees of freedom by ^ in an attempt to adjust for the lack of compound

symmetry. A non-significant unadjusted test implies a non-significant adjusted test

result. For the purposes of exposition, however, we shall complete the analysis.

---

24

The value for the test statistic outlined above is computed below

)

13

)(

81

(

8

68

v

= 273 ^ = 0.6138

A

=

= 0108

.

0

f

=

= 34

M

= 133.23

(

6

)(

273 7)( )

68

2

Q

=

(1-

A

)

M =

131.80

Since

Q

= 131.80 > 2

= 44 the null hypothesis is rejected (

p

= 0.000) and, due

0 05

. ,34

to the lack of compound symmetry, the usual

F

tests will not be valid. If the

calculated value of

Q

exceeds the tabulated chi-square value, then a multivariate

analysis should be used on the data (Milliken & Johnson, 1984).

When the repeated measures are over time, the errors can often be correlated through

an autoregressive error structure. The covariance matrix corresponding to errors with

a first order regressive process is

2

p

-

1

1

L

p

-2

2

1

L

=

2

p

-

1

3

1- 2

L

M

M

M

O

M

p

-1

p

-2

p

-3

1

L

For this data set, the correlations among the residuals from the previous ANOVA are

as follows

---

25

MONTH 2 3 4 5 6 7 8 9

2

1.000

0.597 0.351 0.127 -0.605 -0.254 -0.312 -0.211

3

0.597

1.000

0.669 -0.043 -0.799 -0.109 -0.531 -0.384

4

0.351 0.669

1.000

0.243 -0.726 0.060 -0.545 -0.409

5

0.127 -0.043 0.243

1.000

-0.358 -0.258 -0.529 -0.386

6

-0.605 -0.799 -0.726 -0.358

1.000

0.260 0.742 0.329

7

-0.254 -0.109 0.060 -0.258 0.260

1.000

0.487 0.713

8

-0.312 -0.531 -0.545 -0.529 0.742 0.487

1.000

0.584

9

-0.211 -0.384 -0.409 -0.386 0.329 0.713 0.584

1.000

The correlations are certainly not uniform. For example, the correlation between

Month 3 and Month 6 is -0.799, while the correlation between Month 4 and Month 7

is only 0.060. This disparity is quite large, and not indicative of the variables being

equi-correlated. Also interesting is the way the correlations between Month 9 and the

previous months do not slowly decrease. Since some of the sample residual

correlations fall off with separation in time, an autoregressive error structure is a

possibility.

When the responses of a single subject are measured sequentially, the errors are often

correlated through an autoregressive structure rather than an error structure satisfying

the compound symmetry structure (Lindsey, 1993). Rather than test for an

autoregressive error structure, which can be difficult with only 8 observations per

series of responses, in addition to the fact that many software packages do not supply

a test for this in a repeated measures context, Milliken and Johnson (1984) suggest

analysing the data using multivariate data analysis methods when the number of

subjects exceed the number of repeated measures. The multivariate approach is often

preferable because it is conceptually simpler in that it imposes no restrictive

assumptions on the error covariance matrix ­ in fact it assumes the matrix to be

unstructured, contrast this to the ANOVA models which have several restrictive

assumptions. As was mentioned previously, a non-significant unadjusted test implies

a non-significant adjusted test result. If the experimenter was following the strategy

outlined above they could comfortably stop now, since the

p

value for group × month

is not significant. For the purposes of comparison, however, the multivariate approach

will be applied in the next chapter.

---

26

5. Multivariate Approach to Repeated Measures

Multivariate analysis of variance (MANOVA) is an alternative to repeated measures

ANOVA in which responses to the levels of the within-subjects variables are regarded

as separate variables. In the study being undertaken here, the MANOVA approach

regards the observations at the

p

= 8 different time points as separate variables. Since

the multivariate approach has less restrictive assumptions than the univariate

approach, i.e. MANOVA imposes no restrictive assumptions on the form of the error

covariance matrix, in fact it assumes that this matrix is unstructured. Recall from

Chapter 4 the complicating feature of the univariate analysis is the fact that the

repeated measures are correlated, violating the usual ANOVA assumption of

independence. The choice between the multivariate approach and univariate repeated

measures ANOVA depends on sample size, power, and whether the statistical

assumptions of the univariate ANOVA are met.

The MANOVA model for comparing

r

population mean vectors is

X

= µ + +

e

ji

i

ji

where

i

= 1, 2, ...

nj

and

j

= 1, 2, ...

r

µ is the overall population mean vector

j

is the

j

th treatment effect

and

e

ji

are the

Np

(

0

, ) residuals.

Alternatively, the MANOVA model can be expressed in matrix form

E

(

X

) = ^

A

Where

X

is a

p

×

N

matrix of observations

A

is an

r

×

N

matrix which indicates group membership

---

27

and

^ is a

p

×

r

matrix of parameters estimated from the data.

Since there are only two groups being compared and

p

> 1, Johnson and Wichern

(1982) advise using a multivariate analogue of the two-sample

t

-test called

Hotelling′s

T

2. Hotelling′s generalised Student

T

2 is uniformly most powerful as well

as being robust to violations of the assumptions of multivariate normality and

homogeneity of variance-covariance matrices (Ito, 1980). The test statistic can be

defined as follows.

If

X

1,1,

X

1,2,

X

1,3, ...

X

1,

n1

is a random sample of size

n

1 from an

Np

(µ1, ) and

X

2,1,

X

2,2,

X

2,3, ...

X

2,

n2

is a random sample of size

n

2 from an

Np

(µ2, )

Then

-1

2

1

1

T

= [

X

-

X

- (µ - µ )]′ +

S

[

X

-

X

- (µ - µ )]

1

2

1

2

1

2

1

2

n

n

1

2

is distributed as

(

n

+

n

- )

2

p

1

2

F

p

,

n

+

n

-

p

1

-

1

2

(

n

+

n

-

p

- )

1

1

2

This is the test statistic which will be used in the profile analysis of the data in the

next section. The variance-covariance matrix

S

used in the above equation is given by

1

N

S

=

(

X

-

X

)(

X

-

X

)′.

N

-1

i

i

i

1

=

---

28

5.1 Profile Analysis

Profile analysis is an application of multivariate analysis of variance in which several

dependent variables are measured on the same scale. Profile analysis is useful where

subjects are measured repeatedly on the same dependent variable, and can serve as an

alternative to the univariate repeated measures ANOVA carried out in Chapter 4.

The kinds of questions profile analysis answers depends on the kinds of research

questions asked. For this study the main questions are

(1) Do different groups have parallel profiles? This is known as the test for

parallelism and is the primary question answered by profile analysis (Timm,

1980). For this study, the hypothesis of parallel profiles translates into asking if

the control group and the vaccine group have the same pattern of means on the

various months. In other words, does the vaccine lead to the same pattern of mean

bug production over the course of the experiment as the control group? This is

equivalent to the interaction hypothesis in the repeated measures ANOVA carried

out in Chapter 4

(2) Regardless of whether or not the two groups produce parallel profiles, does one

group, on average, produce higher numbers of cells than the other, ignoring time?

This is called the levels hypothesis in the literature (Johnson & Wichern, 1982;

Tabachnick & Fidell, 1989). It addresses the same question as the between

subjects main effect, `group′, in repeated measures ANOVA.

(3) The third question addressed by profile analysis concerns the similarity of

response to all months, independent of groups. Does each month elicit the same

response? This question is relevant only if the profiles are parallel. If the profiles

are not parallel, then at least one of them is not flat. The flatness test evaluates

whether mean cell production changes over the period of testing. The flatness test

evaluates the same hypothesis as the within subjects main effect, `month′, in

repeated measures ANOVA.

---

29

As with the univariate approach, profile analysis has certain assumptions, all of

which are common to the various forms of multivariate analysis of variance.

5.2 Assumptions of Profile Analysis

Although using a multivariate approach circumvents the problems caused by the

correlated observations, there are some assumptions that should be checked.

The primary requirement of multivariate analyses is that the number of experimental

units in the smallest group, if the groups have unequal numbers of subjects,

exceeds

the number of repeated measurements,

p

. For this study there are 18 cows in the

smallest group and

p

= 8 months, thus this most important of assumptions is safely

met. If this assumption were not met, then the residual matrix would be singular and

the analysis could not be done. In the choice between univariate repeated measures

and the multivariate approach this is often the deciding factor (Tabachnick & Fidell,

1989; Lindsey, 1993). As with the univariate approach, unequal sample sizes provide

no special difficulty.

The second assumption requires that the dependent variables, months in this case,

must have been subject to the same scaling metric. In situations where profile analysis

is used as an alternative to repeated measures designs this is automatically met since

the dependent variables, months, measure the same characteristic but at different

times.

As with the univariate approach, the results of profile analysis generalize only to the

populations from which the subjects are randomly sampled. Since the cows in this

study were randomly selected, this assumption is met.

Profile analysis is as robust to violation of multivariate normality as other forms of

MANOVA. Unless there are under 20 subjects in the smallest group and highly

unequal numbers of subjects in the groups, the assumption of multivariate normality

---

30

is not likely to be violated (Tabachnick & Fidell, 1989). Although there are less than

20 cows in the control group, the sample sizes are not highly disparate.

Multivariate methods are extremely sensitive to outliers. Since this data set contains

numerous outliers, it is important that a test for multivariate normality be carried out.

Vonesh and Chinchilli (1997) recommend carrying out a test for multivariate

normality on the residuals rather than the original observations, while Johnson and

Wichern (1982) advise testing the dependent variables.

Homogeneity of the variance-covariance matrices is important, especially so if the

sample sizes are notably discrepant. Since the sample sizes in this study are not too

disparate, this assumption need not be strictly observed. Univariate homogeneity of

variance is also assumed but the robustness of univariate ANOVA generalizes to

profile analysis.

First, the test of multivariate normality will be carried out.

5.3 Testing For Multivariate Normality

The Multivariate normal density is a generalisation of the univariate normal density.

The univariate normal density with mean µ and variance 2 has the density function

1

-[( -µ ) / ]2 / 2

f

(

x

) =

x

e

2

2

Which leads to the joint density for a multivariate normal population with mean

vector µ and covariance matrix . Assuming that the

p

× 1 vectors

X

1,

X

2, ...

X

N

represent a random sample from a multivariate normal distribution, are mutually

independent, and each has a distribution

Np

(µ, ).

The joint density for all the observations is the product of the marginal normal

densities. Therefore, the joint density of the vectors

X

1,

X

2, ...

X

N

is

---

31

N

1

-

(

x

µ )′ 1

(

x

µ ) / 2

-

-

-

i

i

e

p

/ 2

1/ 2

i

=1 (2 )

| |

N

-

1

-

(

x

-

1

µ )′ (

x

-µ ) / 2

=

i

i

i

=

e

1

Np

/ 2

N

/ 2

(2 )

| |

As previously mentioned, testing for multivariate normality can be carried out on the

original observations or the residuals. Johnson & Wichern (1982) test for multivariate

normality using the original observations, and the method they use will be applied to

the transformed (log

e

) observations.

For each of the

N

= 41 cows, the squared generalised distance

d

2 = (

X

-

X

)′ -1

S

(

X

-

X

),

i

= ,

1 ,

2 ,...,

3

N

i

i

i

is calculated, while the usual unbiased estimator of is

S

. When

p

variables are

observed on each experimental subject, the variation is described by the sample

variance-covariance matrix

s

s

s

1

,

1

,

1 2

L

,

1

p

s

s

s

2 1

,

2,2

L

2,

p

S

=

M

M

O

M

s

s

s

p

1

,

p

,2

L

p

,

p

given by

1

N

S

=

(

X

-

X

)(

X

-

X

)′

N

-1

i

i

i

1

=

---

32

which estimates the population variance-covariance matrix

1

,

1

,

1 2

L

,

1

p

2 1

,

2,2

L

2,

p

=

M

M

O

M

p

1

,

p

,2

L

p

,

p

This definition of the sample variance-covariance matrix is commonly used in many

multivariate test statistics (Johnson & Wichern, 1982). It contains

p

variances and

1

p

(

p

- )

1 potentially different covariances.

2

When the parent population is multivariate normal and both

N

and

N - p

are greater

than 25, each of the

d

2 should behave like a 2 random variable. The values of

d

2 can

2

1

2

1

1

be plotted against ((

i

- ) /

N

) where ((

i

- ) /

N

) is the

i

(

100 - ) /

N

p

2

p

2

2

percentile of the 2 distribution with

p

degrees of freedom. To construct the plot the

values of

d

2 should be ordered from smallest to largest as

2

2

1

2

2

2

2

d

d

d

d

.

(

d

, ((

i

- ) /

N

are plotted. The

)

1

(

(2)

(3)

L

Next, the pairs

))

(

N

)

(

i

)

p

2

chi-squared plot should resemble a straight line if the assumption of multivariate

normality holds. Another way of evaluating

p

-variate normality is if roughly half of

the

d

2 are less than or equal to 2

( .

0 )

5 .

p

Another benefit of the chi-square plot is the ease of identification of outliers. The plot

does not show the presence of any outliers. This is important because multivariate

analysis is extremely sensitive to outliers. Figure 7 shows the chi-square plot for this

data set.

---

33

Figure 7.

20

q8 10

h

iS

C

0

2

7

12

17

DsqSorted

Since the Figure 7 shows little deviation from linearity, the assumption of multivariate

normality is acceptable for the transformed data. The next assumption to be tested is

that of the equality of the variance-covariance matrices.

5.4 Testing For The Equality of Covariance Matrices

The MANOVA approach assumes that the observation vectors for each individual

arise from a multivariate normal distributions, and that the distributions for each

group have the same covariance matrix. The latter assumption is an extension of the

equal variance assumption in univariate ANOVA. The error rates of tests are less

affected if the assumption of equal covariance matrices is false when the sample sizes

are approximately equal. Departures from variance homogeneity have more serious

effects on the Type 1 error rate than departures from multivariate normality (Vonesh

& Chinchilli, 1997). Fortunately, a violation of variance homogeneity has minimal

impact if the groups are of approximately equal size (if the largest group divided by

the smallest group is less than 1.5 (Tabachnick & Fidell, 1989)). A procedure to test

for the equality of variance-covariance matrices is given in Hand and Crowder (1996).

It is less complicated than the more common Box′s M test, and will be used here.

---

34

To test the equality of variance-covariance matrices

H

0 : control = treatment

H

1 : control treatment

Let

S

p

be the unbiased estimator of

p

, with

S

control and

S

treatment estimating control and

treatment, respectively. The likelihood ratio test statistic is

r

r

M

= (

n

- )

1 ln

S

-

(

n

)

1 ln

S

j

p

-

j

j

j

=1

j

=1

This leads to

Mk

following a 2 distribution with (

r

­1)

p

(

p

+ 1)/2 degrees of freedom,

where

k

is defined as

2

p

2 + 3

p

-1

r

1

1

k

= 1-

(

6

p

+ )(

1

r

-

r

)

1

j

=1 (

n

1

j

)-

-

(

n

-1

j

)

j

=1

where

r

= the number of groups in study, 2 in this case. For the data set in question,

the above calculations result in

Mk

= (70.192)(0.78) = 54.7. Since

Mk

is distributed

as 2

with 36 degrees of freedom, the

p

value of the test statistic = 0.03. This test is

extremely sensitive, and provided the sample sizes are not too discrepant, the

hypothesis of equal variance-covariance matrices should only be rejected at the 0.01

level of significance (Hand & Crowder, 1996). Therefore the null hypothesis is not

rejected in this case.

---

35

5.5 Hypothesis Tests for Profile Analysis

Since none of the assumptions appear to have been violated, the profile analysis can

proceed. Both Johnson and Wichern (1982) and Tabachnick and Fidell (1989)

recommend profile analysis as an alternative to repeated measures ANOVA.

The hypothesis tests in profile analysis are the same as those tested in the univariate

ANOVA carried out in Chapter 4, although the terminology is slightly different.

The first hypothesis test is called the test for parallelism in the profile analysis

literature, and as was mentioned earlier, is essentially a test for the interaction

between the groups and the repeated-measures factor. It tests whether the mean

responses along the levels of a repeated-measures factor (months in this case) are

similar for both groups. Parallelism exists when the interaction is not significant. The

mean cell production over time for both groups is shown in Figure 8.

Figure 8.

ln(CELLS)

7

6

5

4

2

3

4

5

6

7

8

9

MONTH

Control

- - - - - Treatment

---

36

Formally, the test for parallelism (interaction between time and group) can be written

as

µ

- µ

0

control

,2

treatment

,2

µ

- µ

0

control

,3

treatment

,3

H

µ

- µ

0

0:

=

control

,4

treatment

,4

M

M

M

µ

- µ

0

control

,

p

treatment

,

p

are true simultaneously,

against

H

1 : at least one µ

control,k

- µ

treatment,k

0, for all

k

from

k

= 2, 3, ...,

p.

Or alternatively,

H

0 : µ

C

= µ

C

1

2

H

1 : µ

C

µ

C

1

2

Where

C

is the (

p

­1) ×

p

contrast matrix

-

0

0

0

0

0

0

1

1

-

0

0

0

0

0

0

1

1

-

0

0

0

0

0

0

1

1

C

=

-

0

0

0

0

0

0

1

1

-

0

0

0

0

0

0

1

1

-

0

0

0

0

0

0

1

1

-

0

0

0

0

0

0

1

1

and µ and µ represent the

p

× 1 mean population vectors of cell production for the

1

2

control and treatment groups respectively.

---

37

The null hypothesis of parallel profiles is rejected at level if

-1

1

1

2

2

T

= (

X

-

X

)′ ′

C

+

′

CSC

(

C X

-

X

) >

c

1

2

1

2

n

n

1

2

n

+

n

-

p

-

2

(

2)(

)

1

where

1

2

c

=

F

()

p

- ,

1

+ -

1

n n

2

p

n

+

n

-

p

1

2

and

X

and

X

are estimates of µ and µ . As usual the population variance-

1

2

1

2

1

N

covariance matrix is estimated by

S

=

(

X

-

X

)(

X

-

X

)′

N

-1

i

i

i

1

=

The transformed data has mean vectors

X

for the control group and

X

for the

1

2

treatment group

373

.

5

633

.

5

127

.

6

088

.

6

803

.

5

048

.

6

1

656

.

5

2

1

n

974

.

5

X

= 1

n

X

=

X

=

X

=

2

1

i

n

945

.

5

i

n

028

.

6

1

i

=1

2

i

=1

978

.

5

087

.

6

342

.

6

129

.

6

373

.

6

118

.

6

---

38

then (

X

-

X

) is

1

2

- 260

.

0

039

.

0

- 244

.

0

- 317

.

0

X

-

X

=

1

2

-

082

.

0

-

108

.

0

212

.

0

255

.

0

For this data set, the calculated value of

T

2 = 10.098, while

c

2 = 19.27 with = 0.05.

The

p

-value of

T

2 can be calculated, with some rearrangement of the above equation,

as

p

= 0.32. Given this result the null hypothesis cannot be rejected, i.e. mean changes

in cell production over time do not depend on group membership. This result is in

accord with the results of the interaction hypothesis tested in the univariate ANOVA.

The next hypothesis test in this profile analysis tests the same hypothesis as the

`Groups′ test in the univariate ANOVA of Chapter 4. It tests whether the group means

differ when the levels of time are collapsed. The mean cell production, ignoring time,

for both groups is shown in Figure 9.

Figure 9.

Boxplots of ln(CELLS) by Group

(means are indicated by solid circles)

8

ln(CELLS)

7

6

5

4

3

2

Control

Treatment

Group

---

39

Thus the levels hypothesis asks whether the means of the two groups are the same,

ignoring times. Formally, the levels hypothesis is

H

µ

= µ

0:

control

treatment

H

µ

µ

1:

control

treatment

This is similar to the univariate two-sample

t

-test. The hypothesis of equal group

means is rejected at the level of significance if

2

1 X

-

2

′(

X

)

1

2

T

=

>

F

()

,

1

n

+

n

-2

1

1

1

2

+

′

1 S

1

pooled

n

n

1

2

For this data set, the value of

T

2 = 0.172, while the = 0.05 level results in a value of

4.08 for the right-hand side of the equation and the

p

-value of this test statistic is 0.68.

Therefore the hypothesis of equal group means cannot be rejected, which is in

agreement with the univariate ANOVA test of group equality.

The final hypothesis tests the main effect for the repeated-measures factor time.

It is a test of whether the whole sample (ignoring group membership) is flat or has a

"profile" (i.e., shows reliable differences across the levels of the repeated-measures

factor ­ here time). This is analogous to the univariate ANOVA test of time effect in

Chapter 3, and in the multivariate literature it is known as the test of flatness. Johnson

and Wichern (1982) liken it to a multivariate generalisation of a one-sample

t

-test.

Figure 10 shows the mean ln(cell) counts for each month, ignoring group

membership. The plot shows that, in general, the mean ln(cell) count increases over

time.

---

40

Figure 10.

6.3

ln(CELLS)

6.2

6.1

6.0

5.9

5.8

5.7

5.6

5.5

2

3

4

5

6

7

8

9

MONTH

Formally, the hypothesis can be defined as

H

µ = µ = µ =

= µ

0:

2

3

4

L

p

H

µ

1: at least two of the

′s are different

k

Or, in terms of the contrast matrix

C

H

0 :

C

µ =

0

H

1 :

C

µ

0

The null hypothesis is rejected at the level of significance if

(

n

+

n

)

X

′ ′

-

C CSC CX

>

F

1

2

[

]′ 1

( )

p

- ,

1

+

-

1

n

n

2

p

---

41

The mean vector of transformed cell counts, ignoring group membership, is

519

.

5

106

.

6

940

.

5

N

834

.

5

X

= 1

X

=

i

N

992

.

5

i

=1

039

.

6

223

.

6

230

.

6

which leads to a test statistic of 32.68 and a critical value of 2.33. Since the

calculated value for the left-hand side of the above equation, 32.68, is larger than the

critical value of the right hand side, 2.33, the hypothesis of µ = µ = µ =

= µ

2

3

4

L

p

can be rejected at the = 0.05 level of significance. The

p

-value of this hypothesis is

0.0000.

In summary, the profile analysis confirmed the results produced by the univariate

ANOVA in Chapter 4. While profile analysis is applicable in many experimental

situations, it has several limitations if an experimenter wants to analyse and fit growth

curves to the average growth of a population over time (Morrison, 1976). The

generalised multivariate analysis of variance is a simple extension of the multivariate

approach used in this chapter, and is the subject of the next chapter.

---

42

6. The Generalised Multivariate Analysis of

Variance

At first it might appear appropriate that linear regression could be used to model the

mean responses for the two groups. The problem, however, is that concerning the

assumptions of linear regression, i.e. the random variables are assumed to be

independent of one another. This is clearly not the case here, as cell counts on the

same cow taken at different times are highly correlated. For example, if the cell count

of a particularly healthy cow in the control sample was low at month 3 relative to the

population cell count at month 3, it would tend to be low at month 4. Thus the

observations on any given cow are not independent across time, so a model other than

the linear regression model may be more appropriate.

6.1 Growth Curves

The Generalised Multivariate Analysis of Variance model (GMANOVA) was

developed by Potthoff and Roy (1964) to fit polynomials or other functions linear in

their parameters to time series of the sizes and weights of organisms. The assumptions

of the GMANOVA are the same as that of profile analysis so the highly correlated

nature of the observations are not problematic.

For the purposes of this study, the model proposed by Potthoff and Roy (1964) and

later by Kshirsagar and Smith (1995) will be used to find the growth curves for both

the control sample and the vaccinated sample. It can also be used to test if the growth

curves are equal for the two groups.

---

43

As was mentioned in the previous chapter, the GMANOVA can be used to

circumvent the problem of correlated observations. Stated briefly, the Generalised

Multivariate Analysis of Variance has the form

E

(X) = B

$

A

where

X

is a

p

×

N

matrix corresponding to

p

= 8

observations on

N

= 41 cows. The

matrix

B

is a

p

×

q

matrix where the rows correspond to the

p

equi-distant time points

and the

q

columns represent the degree of polynomial fitted to the data. The matrix

A

is an

r

×

N

group indicator matrix. The rows of

A

represent the group an observation

belongs to, while each column of

A

represents a subject in the study.

Symbolically,

B

can be represented as

0

1

q

-1

t

t

t

L

1

1

1

0

1

q

-1

t

t

t

L

B

= 1

2

2

L L O L

0

1

q

-1

t

t

t

p

p

L

p

While

A

, the group membership matrix is

1 1 1

0 0 0 0

L

A

=

0 0 0

1 1 1 1

L

Where the rows of

A

represent group membership, here control or treatment, and the

columns represent the

N

= 41 individual cows.

The model that is being developed can be described in more formal terms as follows:

A polynomial regression of the form

E

(

xt

) =

j

0

t

0 +

j

1

t

1 + ... +

jq

-1

tq

-1

---

44

where

t

=

t

1

, t

2,

...

,

tp

;

p

>

q

;

j

= 1, 2,...

r

The notation assumes that there are

r

different groups or treatments and a single

growth variable

x

is measured at

p

time points

t

1

, t

2

, ... tp

on

nj

cows chosen at random

from the

j

-th group (

j

= 1, 2,...,

r

).

For the data being analysed in this thesis, the above notation translates into the

following

expected cell count at time

t

=

E

(

xt

) =

j

0

t

0 +

j

1

t

1

where

q

= 2 and

j

indicates membership of the group 1 (control) or group 2

(treatment) and

r

, which is the number of groups, is 2. Hence

j

= 1 or 2.

Finally, the number of equidistant time points

t

= 2, 3, 4, 5, 6, 7, 8, 9, which are the

months after vaccination the cell counts were taken at, equals 8, which means

p

= 8.

As in the general linear model case, the GMANOVA depends on certain assumptions

being met. These assumptions are the same as that of profile analysis. Since the

assumptions were safely met, there is no need to restate them here.

Since the characteristic being measured is the result of a biological process, a tentative

regression polynomial of

expected cell count at time

t

=

E

(

xt

) =

j

0

t

0 +

j

1

t

1

was chosen to represent the growth curves of the two groups. The regression plots in

the beginning of the thesis indicated that the relationship is close to linear. A better fit

may be obtained by choosing a quadratic term as well, but the plots do not indicate

substantial curvature. In the next section a goodness-of-fit test will be carried out to

ascertain the adequacy of the degree of polynomial.

---

45

The first step in obtaining the coefficients for the growth curves is to obtain the matrix

S

, where

S

is given by

S

=

X I

-

A AA

-

(

′(

′) 1)

X

′

where

X

is the matrix of observations and

A

is the group indicator matrix outlined

earlier.

The above equation results in the matrix

4350

.

17 44

.

1153

.

8 48

.

2 63

.

8 57

.

7 17

.

7 08

.

17

44

.

27 97

.

13 61

.

5 48

.

134

.

9 42

.

4 81

.

4 67

.

1153

.

13 61

.

20 54

.

8 39

.

2 93

.

10 86

.

556

.

4 93

.

8 48

.

5 48

.

8 39

.

15 90

.

4 59

.

7 52

.

4 78

.

4 81

.

S

=

2 63

.

134

.

2 93

.

4 59

.

12 15

.

8 37

.

9 48

.

5 99

.

8 57

.

9 42

.

10 86

.

7 52

.

8 37

.

19 87

.

1141

.

12 92

.

7.17

4 8

. 1

55

. 6

4 7

. 8

9 4

. 8 114

. 1 155

. 0

9 4

. 5

7 08

.

4 67

.

4 93

.

4 81

.

5 99

.

12 92

.

9 45

.

15 98

.

and an estimate of is given by

$

= ( ′ -1 )-1( ′ -1 ) ′(

′)-

B S B

B S X A AA

1

5.80479

5.25771

^ =

0.04205

0.13554

The first column of ^

are the control group coefficients attached to

t

0 and

t

1, while

the second column contains the treatment group′s coefficients. This leads to the

equations

y

= 26

.

5

+ 14

.

0

×

time

for the control and

y

= 80

.

5

+ 04

.

0

×

time

for the

treatment group. Figure 11 shows the two regression lines, while figures 12 and 13

show the connected means and the regression lines superimposed.

---

46

Figure 11.

8

ln(CELLS)

7

6

5

4

2

3

4

5

6

7

8

9

MONTH

Control

- - - - - Treatment

Figure 12.

6.5

)

S

L

6.0

L

E

(

C

ln

5.5

2

3

4

5

6

7

8

9

MONTH

× and o Control

+ and Treatment

---

47

Figure 13.

8

7

6

)

LLS

5

E

l

n

(

C

4

3

2

2

3

4

5

6

7

8

9

MONTH

× and o Control

+ and Treatment

While figures 12 and 13 show the same data, the vertical scale, log

e

(cells) has been

changed to illustrate the often deceptive nature of some plots. While Figure 12

provides ample detail, it gives the reader the impression that the differences between

the groups are quite large. Figure 13′s vertical scale has been resized to reflect the

amount of variability in the data. Recall from Figure 5 that the data were quite

variable, with values between 3 and 8 being plotted. When viewed in the context of

the amount of variability in the data, the differences between the two groups are not

as striking. Transforming the regression lines back to the original scale by applying a

transformation of

ex

to the slope-intercept parameters allows a visual comparison on

the original metric (see Figure 14).

---

48

Figure 14.

1000

S

L

L

500

CE

0

2

3

4

5

6

7

8

9

MONTH

Control

- - - - - Treatment

Control regression line = 192.04 + 1.145 ×

month

Treatment regression line = 331.89 + 1.043 ×

month

6.2 Hypothesis Tests For Growth Curves

There are two hypotheses to test in this chapter:

(1) Is the degree of polynomial adequate to fit the model?

(2) Are the growth curves for the two groups equal?

The first hypothesis is important from a model building perspective. If the polynomial

regression that is fit to the growth curves is not adequate, i.e. the degree of polynomial

is too small, then the amount of unexplained error in the model may be too great.

Conversely, in the interests of parsimony, the degree of polynomial should not be too

great.

---

49

Many researchers (Graybill, 1976; Morrison, 1976) feel that a regression of linear,

quadratic or cubic degree on the time variable

t

is adequate for most growth curves

encountered in a practical setting. From looking at the graphs, a growth curve of

degree one was chosen to initiate the analysis.

If the degree of polynomial is found to be inadequate using the goodness-of-fit test, in

the sense that the true model is quadratic or higher order, then a possible strategy is to

continue to increase the order of the polynomial until the goodness-of-fit test does not

lead to a rejection of the null hypothesis. This strategy is dubious however, because

the final analysis of the data would be based on a series of preliminary tests and this

could enlarge the overall significance level (Vonesh & Chinchilli, 1997). They

suggest sticking with a straight profile analysis model if this is the case.

6.3 Testing Polynomial Adequacy

The hypothesis that a polynomial of degree

q

-1 = 2 - 1 is adequate to describe the

model can be stated formally,

H

0: The degree of polynomial is adequate.

H

1: The degree of polynomial is not adequate.

The level of significance chosen for testing this hypothesis is = 0.01.

A matrix

D

of order

p

× (

p

-

q

) needs to be obtained. The matrix

D

needs to satisfy the

condition

D

′

B

=

0

---

50

where

B

is the matrix

1 2

1 3

1 4

1 5

B

=

1 6

1 7

1 8

1 9

where the first column corresponds to

t

0, the second,

t

1. Naturally, the rows of

B

are

the

p

time points, 8 in this study, at which the cell growth is measured at.

The matrix

D

may be derived by choosing any (

p

-

q

), here (8 - 2) = 6 linearly

independent columns of

I

-

B B B

-

( ′ ) 1

B

′

p

where

Ip

is an 8 × 8 identity matrix. For this study the matrix

D

was found to be

-

0.583333

-

0.333333

-

0.250000

-

0.166667

-

0.083333

0.000000

-

-

0.726190

0.333333

-

0.214286

-

0.154762

-

0.095238

0.035714

-

-

0.250000

-

0.821429

0.214286

-

0.142857

-

0.107143

0.071429

-

-

0.166667

-

0.154762

-

0.869048

0.142857

-

0.119048

0.107143

D

=

-

-

0.083333

-

0.095238

-

0.107143

-

0.869048

0.119048

0.142857

-

0.000000

-

0.035714

-

0.071429

-

0.107143

0.821429

0.142857

-

0.023810

0.083333

-

0.035714

-

0.095238

-

0.154762

0.214286

-

0.000000

0.083333

0.166667

-

0.083333

-

0.166667

0.250000

where the columns of

D

are the first six columns of the matrix

I

-

B B B

-

( ′ ) 1

B

′ .

p

---

51

The GMANOVA test for polynomial adequacy

Source

df

SS & SP matrix

-1

H0

r

H

=

D

′

X

[

A

′(

AA

′)

A X

] ′

D

0

Error

N - r

E

=

D

′

SD

0

Total

N

H

+

E

=

D

′

XX

′

D

0

0

where the appropriate test statistic is Wilks′ Lambda given by

E

0

=

0

E

+

H

0

0

Carrying out the above computations results in

H

0,

E

0 and

H

0 +

E

0 as follows

912

.

1

- 604

.

2

- 570

.

0

981

.

0

253

.

0

450

.

0

- 604

.

2

322

.

4

538

.

0

- 984

.

1

- 492

.

0

-

008

.

1

- 570

.

0

538

.

0

243

.

0

- 093

.

0

- 031

.

0

-

013

.

0

H0

=

981

.

0

- 984

.

1

- 093

.

0

046

.

1

252

.

0

561

.

0

253

.

0

- 492

.

0

- 031

.

0

252

.

0

061

.

0

134

.

0

450

.

0

- 008

.

1

- 013

.

0

561

.

0

134

.

0

307

.

0

333

.

12

- 214

.

6

- 437

.

6

- 269

.

3

- 414

.

1

- 800

.

0

- 214

.

6

223

.

10

- 085

.

0

- 344

.

3

- 772

.

1

781

.

0

- 437

.

6

- 085

.

0

245

.

9

583

.

0

- 178

.

1

011

.

1

E

0 =

- 269

.

3

- 344

.

3

583

.

0

478

.

10

125

.

1

- 626

.

1

- 414

.

1

- 772

.

1

- 178

.

1

125

.

1

275

.

8

-

758

.

1

- 800

.

0

781

.

0

011

.

1

- 626

.

1

- 758

.

1

706

.

4

---

52

246

.

14

- 818

.

8

- 008

.

7

- 288

.

2

- 161

.

1

- 350

.

0

- 818

.

8

546

.

14

453

.

0

- 330

.

5

- 263

.

2

- 227

.

0

- 008

.

7

453

.

0

488

.

9

490

.

0

- 209

.

1

998

.

0

H

0 +

E

0 =

- 288

.

2

- 330

.

5

490

.

0

525

.

11

377

.

1

- 065

.

1

- 161

.

1

- 263

.

2

- 209

.

1

377

.

1

336

.

8

-

623

.

1

- 350

.

0

- 227

.

0

998

.

0

- 065

.

1

- 623

.

1

013

.

5

The determinant of

E

0 = 1690.51, while the determinant of

E

0 +

H

0 = 2806.56,

therefore the test statistic for testing if the degree of polynomial is adequate is

51

.

1690

=

= 6023

.

0

0

56

.

2806

To find the critical value for the hypothesis test the following are required

dm

= order of the error matrix

E

0

dE

= degrees of freedom associated with the error matrix

E

0

dH

= degrees of freedom associated with the hypothesis matrix

H

0,

Since the error matrix is a square matrix of size 6 × 6, it has order 6. The degrees of

freedom for the error matrix are given in the table as

N

-

r

, which for this study results

in 41 - 2 = 39. The degrees of freedom associated with the hypothesis matrix are

given in the table as

r

, which in this study is equal to 2.

Therefore,

dm

= 6

dE

= 39 and

dH

= 2

Kshirsagar and Smith (1995) supply the general rule for rejecting a null hypothesis

based on Wilks′ using a significance level of .

---

53

If

dH

=2 reject the null hypothesis when

1 -

(

d

+

d

-

d

-

H

E

m

)1

×

F

(; 2

d

,

.

m

(2

d

+

d

-

d

-

H

E

m

)1

dm

Therefore, for the test of degree of polynomial with = 0.01, the above equation

results in

1-

6023

.

0

(2 + 39 - 6 - )1

×

F

(

68

,

12

;

01

.

0

)

6023

.

0

6

Since F(0.01; 12, 68) = 2.5 and the left hand side of the equation is equal to 1.635,

the null hypothesis cannot be rejected (p = 0.1038). As a result the polynomial of

degree 1 is an adequate model to represent the growth curves in this study.

As the null hypothesis that a polynomial regression of the form

E

(

xt

) =

j

0

t

0 +

j

1

t

1

cannot be rejected at the = 0.01 level of significance, the next stage of the study will

test the second hypothesis, namely that of the equality of the growth curves for the

control group and the group receiving the vaccination.

6.4 Testing For The Equality of Growth Curves

Formally, the test of the equality of the growth curves for the control group and the

group receiving the vaccination is

H

0: 1

-

2

= 0

---

54

H

1: 1

-

2

0

where 1

is a 2 × 1 vector of coefficients for the control group′s growth curve and 2

is a 2 × 1 vector of coefficients for the vaccination group′s growth curve. In this case,

^

^

^ =

^ =

2,0

1

,10

2

^

^

1

,

1

2 1,

Therefore the hypothesis test is testing whether 1

-

2

= 0

can be represented as

-

,

1 0

2,0

0

=

-

1

,

1

2 1

,

0

To test the general linear hypothesis of the equality of growth curves, the hypothesis

test needs to put into a form

H

0:

L

M = 0

H

1:

L

M

0

where

L

is a matrix of order

l

×

q

and

M

is a matrix of order

r

×

m.

For the hypothesis

under question the choice of matrix

L

is

1 0

L

=

0 1

The columns of

L

correspond to the coefficients attached to

t

0 and

t

1, Therefore

L

is

constructed by putting a 1 in the appropriate column of

L

and 0 elsewhere, depending

on which coefficient or coefficients are of interest.

The rows of

M

correspond to the

r

groups in the study. Because there are two groups

in this study,

M

has the form

---

55

1

M

= -

1

where the columns of

M

represent the contrasts one wishes to compare. Therefore the

equation

L

M

has the form

1 0

1

,

1 0

2,0

-

,

1 0

2,0

×

×

=

0 1

1

1

,

1

2 1

,

-

-

1

,

1

2 1

,

where the matrix

L

M

is estimated by

^

- ^

,10

2,0

^

- ^

1

,

1

2 1

,

which is equivalent to the estimator of the matrix

1 -

2

as given earlier. The matrix

R

needs to be derived using the equation

R

=

AA

-1

I

+

AX S

-1 -

S

-1

B B S

-1

B

-1

B S

-1

XA AA

-

(

′) [

′(

( ′

)

′

)

′(

′) 1]

The MANOVA test for the equality of growth curves

Source

degrees of freedom

SS and SP matrix, order

l

H1

m

H

1

Error

N - r -

(

p - q

)

E

1

Total

N - r -

(

p - q

)

+ m

H

1

+ E

1

Where the error matrix

E

1 is of the form

---

56

E

=

L

(

B

′ -1

S B

) 1

-

L

′

1

and

H

1 is of the form

^

1

-

^

H

= (

L

M

)(

M

′

RM

) (

L

M

)′

1

Carrying out the above equations produces the matrices

R

,

E

1,

H

1, and

H

1

+ E

1 as

follows

0.0086996

0.0839810

R

=

0.0486768

0.0086996

.

23 7637

- .

2 9487

E

=

1

- .

2 9487

.

0 5076

.

2 59674

- .

0

44373

H

=

1

- .

0

44373

.

0

07583

.

26 3604

- .

3 3924

H

+

E

=

1

1

- .

3 3924

.

0 5834

The test statistic for the hypothesis

L

M = 0

is Wilks′ Lambda

E

1

=

1

E

+

H

1

1

For the current hypothesis, this results in

= 3.36658 = 0.870

1

3.86979

---

57

Kshirsagar and Smith (1995) give the procedure for ascertaining the critical values to

test the above.

If

dH

= 1 then reject the null hypothesis when

1-

d

+

d

-

d

H

E

m

×

>

F

(;

d

,

d

+

d

-

d

)

m

H

E

m

dm

Since

dH

, the degrees of freedom for the hypothesis matrix

H

0, is equal to 1, the above

critical value takes the form, for = 0.05

1- 87

.

0

1+ 33 - 2

×

>

F

(

;

05

.

0

,

2

)

32

87

.

0

2

Therefore reject

H

0 if 2.39 > 3.32

Since 2.39 < 3.33 the null hypothesis cannot be rejected,

p

= 0.11, and conclude that

the growth curve for the control group is not significantly different than that of the

treatment group. No analysis if complete without examining the residuals. The

GMANOVA on the transformed (log

e

) cell counts produced residuals which appeared

quite satisfactory. Figure 16 (a) shows the residuals plotted against the predicted

values. The residuals appear uniform over time, with only one outlier. Figure 15 (b) &

(c) show the normal distribution of the residuals.

Figure 15.

Histogram of RES, with Normal Curve

3

50

2

40

1

y

30

0

S

nc

RE -1

20

r

eque

F

-2

10

-3

-4

0

5.5

6.0

6.5

-4

-3

-2

-1

0

1

2

3

FITS

RES

---

58

(a)

(b)

Normal Probability Plot

.999

.99

.95

.80

y

i

lit

.50

.20

r

obab

P

.05

.01

.001

-3

-2

-1

0

1

2

RES

Average: -0.0386646

Anderson-Darling Normality Test

StDev: 0.744071

A-Squared: 0.525

N: 328

P-Value: 0.179

(c)

---

59

7. Conclusion

In this thesis several approaches to the analysis of repeated measures data have been

illustrated. First, the time-by-time ANOVA, a rather simple approach consisting of

conducting

t

-tests of control versus treatment at each month, was considered. In none

of the 8 months of the experiment did the mean cell production differ between the 2

groups. Although the time-by-time approach is not recommended for reasons given in

Chapter 3, it is valuable from a heuristic perspective.

The repeated measures ANOVA provided a more useful, albeit more complex, model

which allows for the fact that the 8 measurements on each cow are correlated.

Univariate models such as this have the advantage of being more powerful than

corresponding multivariate tests illustrated in Chapters 5 and 6 (Vonesh & Chinchilli,

1997). Unfortunately, the assumptions regarding the univariate tests require that the

correlations between the 8 measurements are all the same. This is often an unrealistic

assumption, as those months that are closer together tend to be more closely

correlated than those months which are further apart. An adjustment developed by

Greenhouse and Geisser (1959) can be made if the correlation matrix does not

conform to the equi-correlation structure. This adjustment is only necessary if the

effects being tested are significant. However, Geisser (1980) recommends using a

multivariate procedure if the number of subjects,

N

, exceeds the number of

measurements made per subject,

p

.

The multivariate alternative to repeated measures ANOVA, profile analysis, is

conceptually simpler in that it imposes no restrictions on the form of the correlation

matrix. The results of the profile analysis were similar to the univariate approach in

that no significant interaction was found between group membership and time.

---

60

The table below summarises the

p

values produced by the univariate and multivariate

approaches

Effect Univariate

Multivariate

Group 0.68 0.68

Month 0.00 0.00

Group × Month

0.23 0.32

Given the results of both the univariate and multivariate models, the vaccine to reduce

the cell production in cows does not appear to be effective.

The generalised multivariate analysis of variance allows separate growth curves to be

developed for each group. This allows the information contained within the data set to

be put in a form which is compact and easily interpretable. Unfortunately, analysts

often avoid the use of the GMANOVA due to their lack of familiarity with the

technique and the lack of readily available software (Kshirsagar & Smith, 1996).

In conclusion, carrying out the profile analysis was easier and quicker than the

univariate ANOVA and the GMANOVA. Provided the main assumption of

N

>

p

holds, profile analysis should be used in subsequent experiments of this kind.

---

61

References

Diggle, P.J., Liang, K.Y., and Zeger, S.L. 1994.

Analysis of Longitudinal Data.

New

York : Oxford University Press.

Geisser, S. 1980. Growth Curve Analysis.

In

Handbook of Statistics, Vol 1.

Krishnaiah, P.R. (Ed.), 89-115. North-Holland, Amsterdam.

Graybill, F. A. 1976

. Theory and Application of the Linear Model

. North Sciute,

Mass. : Duxbury.

Greenhouse, S.W., and Geisser, S. 1959. On methods in the analysis of profile data.

Psychometrika

, 24, 95-112.

Hand, D., and Crowder, M. 1996.

Practical Longitudinal Data Analysis.

London:

Chapman & Hall.

Ito, P.K. 1980. Robustness of ANOVA and MANOVA Test Procedures.

In

Handbook
of Statistics, Vol 1.

Krishnaiah, P.R. (Ed.), 199-225. North-Holland, Amsterdam.

Johnson, R.A., and Wichern, D.W. 1982.

Applied Multivariate Statistical Analysis

.

Englewood Cliffs, NJ: Prentice-Hall.

Kshirsagar, A.M., and Smith, W.B. 1995.

Growth Curves

. New York: Marcel-

Dekker.

Lindsey, J.K. 1993.

Models for Repeated Measurements.

New York: Oxford

University Press.

Milliken, G.A., and Johnson, D.E. 1984.

Analysis of Messy Data

, Vol 1. New York:

Van Nostrand Reinhold.

Morrison, D.F. 1976.

Multivariate Statistical Methods

. New York: McGraw-Hill.

Potthoff, R.F., and Roy, S.N. 1964. A Generalized Multivariate Analysis of Variance

Model Useful Especially for Growth Curve Problems.

Biometrika

51:313-326.

---

62

Tabachnick, B.G., and Fidell, L.S. 1989.

Using Multivariate Statistics.

New York,

NY: Harper-Collins.

Timm, N.H. 1980. Multivariate Analysis of Variance of Repeated Measures

.

In

Handbook of Statistics, Vol 1.

Krishnaiah, P.R. (Ed.), 41-87. North-Holland,

Amsterdam.

Vonesh, E.F., and Chinchilli, V.M. 1997.

Linear and Nonlinear Models for the
Analysis of Repeated Measurements.

New York: Marcel Dekker.

---

63

Appendix A

Cow Month2 Month3 Month4 Month5 Month6 Month7 Month8 Month9 Group

1

273 738 416 840 604 321 398 371 0

2

425 748 567 550 1765

1826

3472

2848

0

3

137 3421

150 174 320 339 219 429 0

4

856 1108

273 540 429 1092

2054

1680

0

5

54 87 270 215 263 288 299 738 0

6

300 450 318 108 835 1048

1098

800 0

7

71 160 356 174 230 319 261 271 0

8

133 216 234 135 646 411 454 561 0

9

242 284 228 548 216 171 295 605 0

10 342 1214

1402

877 642 687 769 1145

0

11 71 221 330 207 156 374 254 322 0

12 1878

788 214 195 158 209 447 355 0

13 356 665 342 404 252 729 859 1114

0

14 539 1025

788 557 689 764 998 601 0

15 109 130 115 118 555 48 1024

105 0

16 184 249 344 286 306 379 599 1177

0

17 69 400 350 201 241 260 474 557 0

18 165 531 403 238 259 263 246 268 0

19 2142

785 2209

1304

400 581 476 383 1

20 646 666 529 538 177 572 290 226 1

21

2501 1136 1000 950 900 1811 1487 1403 1

22 328 519 268 102 279 309 448 501 1

23 200 308 374 158 369 273 451 476 1

24 127 326 160 149 250 154 245 196 1

25 579 662 330 433 400 233 430 392 1

26 182 324 543 429 359 315 413 396 1

27 196 383 305 444 583 465 288 271 1

28 809 1939

1354

183 190 434 426 264 1

29 390 857 739 323 829 930 924 1219

1

30 195 268 488 291 350 603 554 613 1

31 252 675 504 325 685 1149

800 890 1

32 601 703 568 608 648 345 359 254 1

33 11 585 352 367 383 395 461 385 1

34 157 78 64 573 363 175 307 290 1

35 323 277 992 560 601 627 404 496 1

36 109 164 232 339 447 458 517 577 1

37 223 598 521 748 516 662 542 361 1

---

64

38 205 513 638 623 609 425 613 378 1

39 958 308 148 210 229 151 132 768 1

40 127 861 704 599 223 409 547 630 1

41 99 72 101 450 817 724 736 577 1

0 = control

1 = treatment

---