On summary measure analysis of linear trend repeated measures data: performance comparison with two competing methods

Unstructured multivariate approach (UMA)

The UMA handles the measurements in the subject as a vector of multivariate responses and treats time points as levels of a qualitative factor with no order. This approach is restricted in equally spaced time points, balanced data with complete measurements and also assumes the homogeneity of covariance matrices in all the k groups.

where the vector ${\epsilon }_{\mathit{i}\mathit{h}}={\left({\epsilon }_{i\mathit{1}h},\mathrm{...},{\epsilon }_{imh}\right)}^T$ is the vector of error for the ith subject in group h.

The primary hypothesis interest in a profile analysis is the parallelism of the k groups' profiles or no group × time interaction effect. The hypothesis can be constructed as H ₀ : C μ ₁ = ... = C μ _k for an appropriate transformation matrix Cwith rank m-1. If the test of interaction is not significant, the tests of the main effects are not confounded. In order to compute any MANOVA-type test statistics such as Wilk's lambda (Λ), the condition N-k >m-1 is necessary, where N is the total number of subjects. Otherwise, the estimated covariance matrix of the transformed responses would not be non-singular and positive-definite. To test time effect, one can investigate the equality of the m elements of the total mean vector (${\bar{\mathit{\mu }}}_{.}$) using one-sample Hotelling's T ² test on the m-1 differences between adjacent measurements from each subject. Here, the same strategy as the SMA is utilized to test group effect, as it is often more efficient than MANOVA-type tests to compare the groups' mean vectors.

Linear mixed model (LMM)

where ${\mathit{X}}_{\mathit{i}}^{\mathit{T}}$ is an m _i × p fixed-effects design matrix for the ith subject, βis a p × 1 vector of fixed-effects parameters for the population, b _i is a q × 1 vector of random effects for the ith subject, ${\mathit{Z}}_{\mathit{i}}^{\mathit{T}}$ is an m _i × q random-effects design matrix for the ith subject with q ≤ p, and ε _i is an m _i × 1 vector of within-subject errors. The random-effects vectors, b _i , are assumed to be independent and to have a multivariate normal distribution with mean zero and covariance matrix G _i , and the error vectors, ε _i , are assumed to be independent and to have a multivariate normal distribution with mean zero and covariance matrix R _i . In addition, it is also assumed that b _i and ε _i are independent of one another. The LMM defines the covariances of the measurements in the subject by the covariances of the random effects (G _i ) and the covariances of the errors (R _i ). We used the estimators based on the restricted maximum likelihood (REML) method to construct the F statistics of the hypotheses since, in general, it yields less biased estimates of the variance components than those of maximum likelihood (ML) approach and avoids inflating type I error rates [12, 13].

The summary measure approach (SMA)

In this section, we describe how to apply the least squares regression slope and mean of response over time for each subject to test the effects of time, group and group × time interaction in repeated measures studies.

The slope of least squares regression line was applied to summarize the relationship between response and time for each subject or within-subjects effect. If the pattern of individual profiles is linear or at least monotonic, the slopes can appropriately summarize the rate of change of response over time in the subjects. For repeated measures designs, the primary hypothesis is to test whether the pattern of change over time is the same across the k groups or no group × time interaction effect. Under the assumption of no interaction effect, the slopes in the k groups should not be significantly different. For this purpose, once the slopes are obtained for each subject, the ordinary k sample tests such as one-way ANOVA F or Kruskal-Wallis (for k > 2) and Student's t or Wilcoxon-Mann-Whitney (for k = 2) can be employed to assess the equality of the slopes in the groups. If the test of interaction is not significant, one would be interested in assessing the main effects.

The hypothesis of no time (within-subjects) effect states that all the m elements of the total mean vector (${\bar{\mathit{\mu }}}_{.}$) are identical. Under this assumption, the overall mean of the slopes in the population must be zero. To test this hypothesis, one-sample t test can be applied to the sample slopes to assess the departure of mean slopes from zero.

For testing group (between-subjects) effect, the mean of measurements over time for each subject is used as a summary measure. By analogy with the interaction effect case, the ordinary k sample tests are applied, but this time, to assess the equality of the individual means in the groups.

Permutation procedure can also be employed to assess the interaction and group effects where the constructive assumptions of the standard tests are not held or cannot be reasonably checked due to small sample sizes in the groups.

Simulation study

where Y _ij is the jth measurement from ith subject and X _i is a grouping variable with the values 0 and 1, for i = 1,...,N and j = 1,...,m _i .

Linear trend mixed model data was generated based on the model (3) with the same measuring time points t _ij = t _j = 2j for all the subjects, m _i = m = 5, 10, and 20 measurements and β ₀ = 2, in which the random effects, b _0i , were assumed to be independently normally distributed with mean zero and standard deviation 0.25.

Since hypothesis testing effects related to within-subjects effect is highly dependent on the number of measurements, the values of β ₂ and β ₃ are adjusted with respect to the m values. Different combinations of β ₁ , β ₂ and β ₃ were constructed to compute the empirical type I error rates and powers for testing the three effects.

We considered the following three covariance structures for errors to generate artificial data and fit the LMMs:

• Simple or independent (IND): R _i = σ ² I, where Iis an m × m identity matrix.

• First-order autoregressive (AR1) with ρ = 0.7: R _i = σ ² H, where H= [h _jj' ] is an m × m matrix with h _jj' = ρ ^|j-j'|for all j and j'.

• Unstructured (UNS): R _i = [r _jj' ] is an m × m covariance matrix with arbitrary structure.

For simplicity, we defined the true structures as those which were used to generate data and the working structures as those which were used to fit the model. In all the cases, it was assumed that the errors were normally distributed with zero mean and in the cases of IND and AR1, the error variances were fixed over time and equal to σ ² = 0.5.

1000 sample data sets were generated for n ₁ = n ₂ = n ₃ = 5, 10, 30 and 50 subjects under various choices of the above circumstances.

We have used free statistical software environment R to generate the artificial datasets and fit all of the approaches presented in the method section.

Simulation results

Illustrative examples

Example 1: Pituitary-pteryomaxillary distance data

The first example is a small data set on a facial distance previously published by Potthoff and Roy [14] conducted at the University of North Carolina Dental School. The distance (mm) from the centre of the pituitary gland to the pteryomaxillary fissure was measured at age 8, 10, 12, and 14 in two groups of children (11 girls and 16 boys). The data set has also been analyzed by several analytic methods [12, 15].

Figure 1 displays the mean profiles in boys and girls and indicates a departure from the parallelism hypothesis. In general, boys tend to have larger pituitary-pteryomaxillary distances and a faster growth rate than girls. In addition, the distances increase over age points in both groups of children.

Given model (3), three models were fitted with IND, AR1, and UNS covariance structures. Random intercept and slope models with the three covariance structures were also employed. Table 4 reports the results for the six LMMs, UMA and SMA, as well as Akaike's information criterion (AIC) and Bayesian information criterion (BIC) as two model selection indices for the LMMs. The model with the smallest criterion provides the best fit to data. Based on both foregoing criteria, model 1 with IND covariance structure was preferred. It is worth mentioning that a random intercept model with IND covariance structure for errors yields a compound symmetry covariance structure between the responses.

			Effect
Method	no	Structure	Group	Time	Group × Time	AIC	BIC
LMM with Random intercept	1	IND	0.005	< 0.001	0.014	445.76	461.62
	2	AR1	0.006	< 0.001	0.012	447.71	466.22
	3	UNS	0.007	< 0.001	0.009	450.17	481.90
LMM with Random intercept and slope	4	IND	0.009	< 0.001	0.026	448.58	469.74
	5	AR1	0.015	< 0.001	0.021	446.81	470.61
	6	UNS	0.010	< 0.001	0.009	454.12	491.14
SMA	7	-	0.005	< 0.001	0.019	-	-
UMA	8	-	0.005	< 0.001	0.070	-	-

All the LMMs and the SMA showed a significant interaction effect at the 5% significance level. These results indicated that the growth pattern in boys was faster than that in girls. Although one could not reject the hypothesis of no interaction effect by the UMA, there was some evidence that the profiles in Figure 1 were not parallel.

All the approaches yielded significant results for the two main effects on the facial growth measurements of children. Based on these results, we accept that boys have larger facial distances than girls and the facial distances increase over age in the two groups of children.

Example 2: Change in lung NO metabolites level data

The second example is an animal experimental study which is about the effects of hypercapnia with or without acidosis on NO production in the isolated ventilated-perfused rabbit lung by assessment of the NO metabolites (nitrite and nitrate) concentration released into the perfusate. The study was conducted at Justus-Liebig-University, Giessen. The NO metabolites concentration (nmol/min) was measured at time point 0, 5, 10, 15, 30, 45, ..., and 180 minutes in three groups of normoxic normocapnia (NX-NC, n = 7), normoxic hypercapnia with acidosis (NX-HCA, n = 4) and normoxic hypercapnia with normal pH level (NX-HCN, n = 6). Since there were some variations between the baseline measurements, values were given as changes from the baseline. There were six samples (lungs) with incomplete measurements.

Figure 2 displays the mean profiles of change in NO metabolites level data over time for the three groups. The mean profiles increase over time points in all of the groups. However, it is not expected that the patterns of change in NO metabolites level and the overall means will differ between the three conditions.

In this data set, the UMA could not be conducted, because the number of measurements (m = 14) was larger than that of the samples with complete measurements (11 lungs) during the period of study. On the other hand, the UMA and repeated measures ANOVA approaches are not able to handle the experimental units with missing observations.

Table 5 displays the results for the LMM with random intercept, random intercept and slope and also the SMA. Note that AIC prefers random intercept and slope model with UNS, IND and AR1 covariance structures, models 6, 4 and 5, respectively, whereas random intercept and slope model with IND and AR1 covariance structures are to be preferred based on BIC, models 4 and 5, respectively. The reason is that a heavier penalty in the calculation of BIC than AIC was imposed when the number of parameters in the model increased. Since there were a limited number of lungs and a large number of measurements, the danger of over-fitting increases. In these cases, it is more reasonable to rely on BIC to select the best parsimonious model. Note that model 6 has larger parameters (d = 112) than model 4 (d = 8) which must be estimated.

			Effect
Method	no	Structure	Group	Time	Group × Time	AIC	BIC
LMM with Random intercept	1	IND	0.264	< 0.001	< 0.001	121.38	148.31
	2	AR1	0.248	< 0.001	0.026	58.77	89.06
	3	UNS	0.175	< 0.001	< 0.001	10.91	344.14
LMM with Random intercept and slope	4	IND	0.212	< 0.001	0.236	4.51	38.17
	5	AR1	0.209	< 0.001	0.234	5.90	42.93
	6	UNS	< 0.001	< 0.001	0.262	-0.57	339.39
SMA	7	-	0.273	< 0.001	0.241	-	-

Based on the results of the LMMs 4 and 5 selected on the basis of BIC, and also the SMA, one can accept that the rates of NO metabolites change in the three groups do not differ. Although this result coincides with that obtained by the most complicated model 6, the unsuitable models 1 and 3 reject the hypothesis of no interaction effect which is not illustrated in Figure 2.

All the LMMs, as well as the SMA, confirmed the effect of time on increasing the mean change over time in all of the groups. Except for the unreasonable model 6, all the models and the SMA confirmed that the mean change profiles for the three groups were the same throughout the time points; therefore, there was no significant group effect.