Repeated-Measures Analysis in the Context of Heteroscedastic Error Terms with Factors Having Both Fixed and Random Levels

The design and analysis of experiments which involve factors each consisting of both fixed and random levels fit into linear mixed models. The assumed linear mixed-model design matrix takes either a full-rank or less-than-full-rank form. The complexity of the data structures of such experiments falls in the model-selection and parameter-estimation process. The fundamental consideration in the estimation process of linear models is the special case in which elements of the error vector are assumed equal and uncorrelated. However, different assumptions on the structure of the variance–covariance matrix of error vector in the estimation of parameters of a linear mixed model may be considered. We conceptualise a repeated-measures design with multiple between-subjects factors, in which each of these factors has both fixed and random levels. We focus on the construction of linear mixed-effects models, the estimation of variance components, and hypothesis testing in which the default covariance structure of homoscedastic error terms is not appropriate. We illustrate the proposed approach using longitudinal data fitted to a three-factor linear mixed-effects model. The novelty of this approach lies in the exploration of the fixed and random levels of the same factor and in the subsequent interaction effects of the fixed levels. In addition, we assess the differences between levels of the same factor and determine the proportion of the total variation accounted for by the random levels of the same factor.


Introduction
Linear mixed-effects models can handle both fixed and random effects simultaneously. Not only are the models convenient for modelling the means of the data, but they are also convenient for modelling covariances [1,2]. Situations arise in fields such as psychology and medicine in which longitudinal or correlated data display notorious heterogeneity of responses to stimuli and treatment [3]. When the classical linear mixed-effects model is inappropriate, alternative approaches, which identify and consider subgroups in fixed and random effects, are necessary. The authors of [4] propose a heterogeneity model as an extension to that proposed by the authors of [5,6] by replacing the normality assumption for random effects and by assuming a non-zero mean vector and a common error vector.
Repeated-measures data are often dependent, a property that does not conform to the generality of a mixed-effects model [1,7,8]. Linear mixed-effects models are some of the most convenient statistical approaches, which account for this dependency [9]. However, setting them up for data analysis requires some care, especially in choosing the most appropriate covariance structure to keep the type I error down [10]. Selecting a suitable model for the covariance is important because the precision of the confidence intervals and the tests of hypotheses concerning model parameters depend on the correct model [11]. According to [12], selecting an appropriate mixed-effects model and construction approach, such as partitioning of fixed and random effects, allows for analysis of variance in correlated data. Numerous diagonal and non-diagonal covariance structures for correlated data, Stats 2022, 5 which cover a range of assumptions about the associations between responses from the same cluster data, are available [13,14]. Modern statistical software, such as SAS Studio in Linux Environment, provides options for selecting candidate covariance structures through the PROC MIXED procedure [15] (see Appendix A). The residual maximum likelihood (REML) is one of the most famous methods for estimating covariance parameters associated with linear mixed models [16,17], among other alternatives [18].
We present the construction and analysis of a three-factor linear mixed model for repeated-measures designs when the between-subjects factors consist of both fixed and random levels, and the structure of the variance-covariance matrix of the error terms is not the identity. We consider a repeated-measures design setting under a linear mixedeffects model, with factors sharing both fixed and random components of the model. For experimental designs consisting of factors with a unique composition of both fixed and random levels, we propose a partitioning approach based on factor levels in model construction, estimation, and hypothesis testing. In such a situation, the fixed levels allow for the comparison of specific levels of interest within the factor, whereas the random levels allow for the assessment of variation within the same factor. We assess the effect of introducing heterogeneity of error terms in the selection of the most appropriate covariance structure, in the assessment of the changes that occur in the estimation, and when drawing inferences. To simplify the proposed approach, we focus on the diagonal covariance structure for the repeated-measures linear mixed-effects model. Section 2 presents the approach to model construction, for a general linear mixedeffects model in a completely randomised design (CRD) and in a repeated-measures design (RMD) under the default assumption of error terms. Section 3 presents results from a numerical example. Discussion of the results is in Section 4, followed by the conclusions in Section 5.

An Illustrative Data Structure
We motivate the approach using data collected from a study that investigated the impact of combining carbon tetrachloride (CCl 4 ) with four levels (0, 1.0, 2.5, and 5.0 mM) and chloroform (CHCl 3 ) with four levels (0, 5, 10, and 25 mM) on the toxicity of cells on in vitro toxicity of isolated hepatocyte suspensions [19]. Four flasks were assigned to each of the 16 treatments. Cell toxicity is measured by the amount of lactic dehydrogenase (LHD) enzyme percentage leakage from each of the 64 flasks after 0.01, 0.25, 0.5, 1, 2, and 3 h since the application of the treatment. For illustration purposes, we consider the between-subjects factor CCl 4 levels 2.5 and 5.0 as fixed (new technology), and levels 0 and 1.0 are considered as existing random levels (old technology). Similarly, we consider the between-subjects factor CHCl 3 levels 10 and 25 as fixed (new technology), and levels 0 and 5 are taken as random levels (old technology). Of interest in the analysis is the percentage leakage observed at times 1, 2, and 3. We demonstrate the model construction procedure under certain assumptions in a completely randomised design (CRD) and in a repeated-measures design (RMD).

Construction of a Linear Mixed-Effects Model in CRD
Consider a three-way treatment structure in a balanced, completely randomised design (CRD) with full interaction of factors A, B, and C, each consisting of f fixed and r random levels. Assume we have f A , f B , and f C fixed levels and r A , r B , and r C random levels of factor A, respectively. We partition the dataset based on the combinations of factor levels and construct a partitioned model in each partition. For example, the FRF partitioned model is built from the f A , r B , and f C levels. Similarly, other possible partitions are FFF, FFR, RFF, RRF, RFR, FRR, and RRR. We illustrate the model construction using the FRF linear mixed-effects model in CRD, having at least one replication per treatment combination and expressed as where y FRF ijkl is the lth observation in the (ijk)th treatment cell of the FRF partition; l = 1, . . . , r h are the replicates (where all r h = r for balanced data); µ FRF is the overall mean; ϕ Ai , ϕ Bj , ϕ Ck are the main effects of the three factors; π 1 , . . . , π t are the interaction effects; . . , c (c = f C + r C )) are unknown parameters corresponding to fixed factor A, random factor B, and fixed factor C, respectively. Defining the random main effect as ϕ R and the random interaction effect as π R in (1), the random effects and the random error term ijkl , are commonly assumed to have zero mean and variance, i.e., ϕ R ∼ N 0, σ 2 ϕ R , π R ∼ N 0, σ 2 π R , and ijkl ∼ N 0, σ 2 . For a balanced data scenario with r replications per cell, for example, the general linear mixed model Equation (1) is normally expressed in matrix form as where y FRF : N × 1 is a vector of response observations in the FRF partition; matrix X FRF : N × p is a known incidence matrix associated with the vector of p fixed-effects β : p × 1 in the model; matrix Z FRF : N × q is a known incidence matrix associated with the vector of q random-effects u : q × 1 in the model; and FRF : N × 1 is a vector of random errors. The usual assumption under this model is that the random effects are u ∼ N(0, G), and the random residuals are FRF ∼ N(0, R), where R = σ 2 N and G is a diagonal matrix of variance components (i.e., different variances, and all zero covariances).
The total variance-covariance is the structured matrix V = ZGZ + R, a structure that guarantees independence and homogeneity of residual errors. This implies that the variance of y is modelled through Z, G, and R. The simple total variance-covariance, V, has a block-diagonal structure given by the matrix

Linear Mixed-Effects Model in RMD
Traditionally, between-subject and within-subject factors in repeated-measures experiments are designated as either fully fixed effects or random effects. Depending on the objectives of the experiment, some factors in linear mixed-effects models may exist with both fixed and random levels [12,20,21]. The same scenario is common with a repeatedmeasures experiment, in which either the between-subjects factor or the repeated measures consist of both fixed and random levels. For instance, for improved results, a researcher may decide to consider additional levels of a between-subjects factor, in addition to the old and existing levels. In that case, the improved analysis needs to consider the new factor levels as fixed levels, and the old and existing levels are considered random. The approach creates an opportunity to compare and evaluate the effectiveness of new factor levels (fixed) against existing ones (random) and/or to compile a combined analysis of both. In addition, random levels allow for the assessment of variation between and within the factor levels in the entire population.
We consider a three-factor repeated-measures experiment with n experimental units (EU) that are randomly assigned to each of the a levels of the between-subjects factor A (with f a fixed levels and r a random levels, a = f a + r a , for example), the between-subjects factor B (with f b fixed levels and r b be random levels, b = f b + r b , for example), and the within-subjects factor C with t measurements (all considered to be fixed in this case) taken on each of the experimental units (EU).
where the subscript p in y p ijkl and p ijkl denotes the partition (p = 1, 2, . . . , 4) and where y p ijkl is the lth measurement (l = 1, 2, . . . , t) of the kth experimental unit (k = 1, 2, . . . n) in the ith level (i = 1, 2, . . . a) of factor A and the jth level (j = 1, 2, . . . b) of factor B in the pth partition. Depending on whether the effects are fixed or random, the model parameters are as defined in (1).

Model Assumptions
In experimental research that involves analysis of variance (ANOVA) as a technique for comparing different treatment means, a set of assumptions that includes the usual normality and homogeneity of variance must be checked before the analysis of data [22]. These traditional normality tests, such as Q-Q plots or the Shapiro-Wilk test, and outlier detection approaches (e.g., box plots), are appropriate for diagnosing violations of assumptions.

Sphericity (Circularity) Assumption
Similar to the homogeneity of variances in a between-subjects analysis of variance, sphericity holds in ANOVA for repeated-measures designs when the variances of the differences among all possible pairs of the within-subject factor group means are equal [23,24].
The assumption is usually unrealistic in repeated-measures designs in which observations are correlated. Univariate tests for within-subjects effects apply when sphericity holds. If the sphericity assumption is violated, several alternatives for adjusting the numerator and denominator degrees of freedom have been suggested [25]. The most popular degrees of freedom corrections are the Greenhouse-Geisser [26] or the Huynh-Feldt [27] adjustments, which are extensions of the Box correction factor [28]. These adjustments provide a more accurate adjusted p-value [29,30]. Failure to address the problem of sphericity when conducting an analysis of variance often leads to inflated F-ratios, type I errors, biased post hoc tests for group mean differences, and inaccurate conclusions [24]. Tests for sphericity are offered in most statistical computing packages (SAS, R, SPSS, etc.).

Compound Symmetry Assumption
An overly restrictive assumption closely related to sphericity is the compound symmetry assumption. It states that there is a constant variance and correlation between observations from the same subject, which is not always realistic in many repeated-measures applications [31]. The Greenhouse-Geisser [26] or the Huynh-Feldt [27] adjustments are used to circumvent restrictive compound symmetry assumptions and to accommodate more general covariance structures for the repeated measures [1]. Compound symmetry implies sphericity, but not vice versa. Univariate tests for within-subject effects apply when compound symmetry holds. If compound symmetry (and hence sphericity) holds, then a split-plot analysis is used as an appropriate approximation to the repeated-measures experiment, since it provides relatively more accurate p-values for testing treatment effects [31]. Compound symmetry is simplified in the Huynh-Feldt [27] condition for sphericity. Hence, both are tested in many software packages by Mauchly's test [32]. However, Mauchly's test for sphericity has been criticised for its over sensitivity and tendency to reject compound symmetry [33].

Estimation Techniques
The general linear mixed-effects model (2) is used to describe data from partitioned repeated measurements, wherein the fixed-effects component, X p β, consists of the design matrix X p , and the fixed-effects coefficients β are as defined in (3). The random-effects component, Z p u, contains the block-diagonal random-effects design matrix Z p with design matrices for the individual subjects (Z i ); u is a vector of random coefficients (the betweensubject variance-covariance components), and p denotes the within-subject errors from the pth partition. The random effects and the residuals follow the distribution u ∼ N(0, G) and p ∼ N(0, R), respectively, where G is a block-diagonal covariance matrix of the random effects and where R is a diagonal covariance matrix with partitions corresponding to each subject (within-subject errors) in the analysis. The covariance matrix for the repeated-measures data is composed of matrices Z p , G, and R, which is a block-diagonal Σ p = Var y p = Z p GZ P + R. The non-singular components G and R are usually estimated by two principal likelihood methods for estimating variance components [34], i.e., the maximum likelihood (ML) method and restricted maximum likelihood (REML). These procedures are available in various mixed-model statistical software, such as SAS (PROC MIXED procedure) and R, with the REML estimates generally preferred unless the data sets are quite large [34].
Assuming that both the random effects and the error terms are normally distributed, the likelihood function for the repeated-measures mixed model is given by [35]: where y p and V = Σ p are as defined in (2). Similarly, a modification of the ML procedure through factorisation of the likelihood function is proposed as an alternative method of estimating covariance parameters, which is the restricted maximum likelihood function [36]: where the available covariance matrix Σ p is used to estimate the fixed-effects parameters, The main challenge in repeated-measures analysis of variance is to determine the adequate correlation structure, because the constant variance assumption for the distribution of the error terms is likely not to be reasonable for the distribution of error terms within subjects. There are various possible choices of covariance structures for repeated measures within each subject depending on the chosen parameterisation for G and R. The choices are usually guided by limitations of the software and by insight from the researcher. The most common covariance structures include variance components, compound symmetry (common covariance plus diagonal), unstructured (general covariance), and autoregressive [37]. With the PROC MIXED statement in SAS, one can specify any repeated-measurements covariance structure for G by using the Random statement and by specifying the form of R with the Repeated statement, in conjunction with the Type option [38]. Excluding the Repeated statement specifies the classical R, which is assumed to be equal to σ 2 I.
There are numerous ways of identifying the most appropriate covariance structure amongst set candidate structures [37]. The most recommended approach is to select the structure that gives the smallest Akaike's Information Criterion (AIC) [39], a statistic that is defined by the model and the maximum likelihood estimates of the parameters from specifying the variance-covariance as where g is the effective number of independently adjusted parameters in the covariance matrix from the pth partition, and L β ,Σ p = log(ML) is the value of the likelihood function evaluated at β ,Σ p . A better model is the one with the smallest AIC value. Different forms of R can be compared for adequacy using the likelihood ratio test statistic [38]. The hypotheses involved are the following: H 0 : R 1 is as adequate as R 2 ; and H 1 : R 1 is not as adequate as R 2 , where R 1 is a special case of R 2 . Suppose R 1 and R 2 have g 1 and g 2 parameters, respectively, with (g 1 < g 2 ), and the test statistic is

Methods of Inference
We present the algorithm for obtaining expected mean squares using the ANOVA approach for the FRF model. Once the coefficients of variance components are established in step (m) of Figure 1, the expected mean squares are found by Step (h) of the Process in Section 2.6.1. Table 1 summarises the variance components obtained for a three-factor repeated-measures design, with two between-subjects factors (A and B) on one within-subject factor (C) using the Process in Section 2.6.1. We use this process to build an FRF model, which assumes factors A and C as fixed, factor B as random, and experimental units (EU) as random. The F, R, F, and Rep denote the effect; a, b, t, and n denote levels of factors A and B, treatment number, and sample size, respectively. The subscripts i, j, l, and k denote the indices of the corresponding levels. Figure 1 illustrates the key steps to be taken when constructing the coefficients of variance components for a partitioned linear mixed-effects model.

Source of Variation
Factor Effect For example, E(MSA), with effect α i , is composed of the variance components θ α , σ 2 αβ , σ 2 γ(αβ) , θ ατ , σ 2 αβτ , and σ 2 as follows: Table 2 displays the ANOVA layout and the expected mean squares for a three-factor repeated-measures design when one of the factors is a within-subject factor. The FRF model is considered for illustration, with experimental units assumed to be random. Table 2. Expected mean squares for a three-factor repeated-measures design.

Hypothesis Testing for Fixed Effects
We are interested in testing the main and interaction effects of the between-subjects factor and of the within-subjects factor in both the partitioned and combined repeatedmeasures linear mixed-effects model. In addition to checking model assumptions, the following hypotheses are of interest for each partitioned model:
The test statistic for H1 is given by F = MSA MSU(A×B) ∼ F a−1,ab(n−1) (α), and the test statistic for H2 is given by F = MSP MSE ∼ F t−1,ab(t−1)(n−1) (α). The interaction effects among the between-subjects and the interaction effects among within-subjects are tested by the MSU(A × B) and the MSE on the denominator, respectively.

Hypothesis Testing for Random Effects
Variance components are estimated by equating the mean square to the expected mean squares derived in Table 2. Where there are no valid F tests, approximate F tests are constructed for the sources of variability in random effects [40]. For random factor B, the hypothesis of interest may be Hypothesis 3 (H3). Random effects (e.g., H 0 : σ 2 β = 0, against H 1 : σ 2 β > 0).

Combined Analysis
The individual partitioned models provide pieces of information which are needed for an integrated analysis. The combined model is built by combining the degrees of freedom and sum of squares associated with each source of variation for each appropriate hypothesis test. For example, the combined effect of the within-subjects factor (Period) in the FC model is obtained from the partitions in which the factor Period is fixed, i.e., the pieces of information are supplied by the partitioned models FFF, FRF, RFF, and RRF. Similarly, the other main and interaction effects for the combined model are obtained by summing up the associated degrees of freedom and sums of squares.

Checking Model Assumptions
A three-factor (CCl 4 , CHCl 3 , and Time) repeated-measures experiment on one of the factors (Time) is proposed. Table 3 shows the multivariate data (wide format) layout for a three-factor repeated-measures experiment.

Normality and Outliers
Q-Q plots are used to check normality and outlier assumptions simultaneously. The FFF, FRF, RFF, and RRF data sets do not show any serious deviations from normality. Figure 2 shows the normal Q-Q plots for the four data subsets. Furthermore, the Q-Q plots do not show any influential point (outlier) that warrants exclusion, since all plots are not very far from the diagonal.

Sphericity and/or Compound Symmetry
In order to test the suitability of using a repeated-measures design in the experiment, the sphericity (or compound symmetry) assumption is tested for each partitioned data subset. Sphericity test results are produced by fitting two models using the PROC MIXED procedure: one specifying the unrestricted covariance structure and the other with a less conservative Huynh-Feldt (H-F) adjustment in the Type option [41]. The null hypothesis for the test is :   Hypothesis 4 (H4). The covariance structure fits the sphericity structure.
The difference (D) between the two −2 log-likelihoods of the two competing models follows a Chi-square distribution with degrees of freedom equal to the difference in the numbers of parameters in them.

Analysis of Results
We use the restricted maximum likelihood estimation (REML) approach to estimate these variance components. Table 4 contains the estimated AIC from each of the covariance structures as well as the number of covariance parameters estimated with a non-zero value. The covariance structure ARH(1) has the minimum AIC values in models FFF, FRF, and RRF, and covariance structure AR(1) has the smallest AIC value in model RFF. The covariance structure ARH(1) is chosen as the most adequate covariance structure for the partitioned models based on AIC.
When the appropriate covariance structure is incorporated into the model, the relationships among the errors and variance components are specified. The PROC MIXED procedure does not require further assumptions on them [42], making it a robust and flexible procedure [15] that produces valid F-tests regardless of whether the sphericity assumption is satisfied or not. Table 5 summarises the PROC MIXED F-tests for the main and interaction effects of the between-subjects and the within-subjects factors for the four partitioned models. The CHCl 3 factor has a significant effect on the leakage percentage (p-value < 0.001) over time in the FFF partition. Based on the fixed-effects (FFF) repeated-measures study, it may be concluded that chloroform (CHCl 3 ) has a significant impact on the amount of lactic dehydrogenase (LHD) enzyme percentage leakage (toxicity of cells) over time, whereas neither carbon tetrachloride (CCl 4 ) in isolation nor the interaction thereof has significant influence. Furthermore, the Time factor plays an important role in determining the amount of leakage as well. However, the interaction of old (random) and new (fixed) levels of the between-subjects factor levels has non-significant effects in FRF, RFF, and RRF models when CHCl 3 fixed levels are involved. Table 6 gives a summary of the estimated covariance parameters and the proportion (in percentage) of variation they contribute in each of the partitioned models. A zero variance for the random levels of chloroform (CHCl 3 ) and very small estimates in other factors are obtained, which result in very low estimates of covariance parameters in the models. The CHCl 3 and CCl 4 random levels and their interaction have a noticeable contribution to the proportion of variation in the amount of lactic dehydrogenase (LHD) enzyme percentage leakage in the RRF partition. Generally, Time has very little interaction effect with the between-subjects factors in FRF and RRF models in determining the proportion of variation of the toxicity of cells.
Since the mixed-model methodology directly computes neither sums of squares nor F-statistics from the ratio of mean squares, and since the R CRAN package lme4 for mixed models currently does not have options for other covariance structures to cater for correlated error variances, generating a combined analysis is not a straightforward exercise. We scrap the targeted data subset based on the effects of interest before fitting the models. Analogous to the partitioned analyses, a comparison of model fit via the AIC approach is conducted. For convenience purposes, let the factors CCl 4 , CHCl 3 , and Time be labelled as factor A, B, and C, respectively. The PROC MIXED procedure is used to fit the repeatedmeasures linear mixed models for intended narrow inference spaces [12,43]. Table 7 shows the Type III tests for the combined models.
Of the possible candidate covariance structures (CS, CSH, AR(1), and ARH (1)), structure ARH(1) is selected as the most appropriate covariance structure for the combined fixed-effects model FA, and AR(1) is appropriate for FB and FA × FB. The factors B (CHCl 3 ) and C (Time) have significant effects (p-value < 0.05) in the combined models. The broad inference scope results for the combined models (assuming random factor A or B effects) are similarly analysed. Table 7. Fixed-effects F-tests for the combined models in narrow inference scope.

Discussion
Based on the illustrative example results, the approach manages to isolate the effects of new and old factor levels over time. The combined analysis confirmed the results of the partitioned analysis on the percentage leakage in cells. The proposed approach conforms to the model construction and the analysis procedures in a repeated-measures design. It can be used as a planning tool in which factor combination and time are of interest in designing experiments that involve repeated measures. In such experiments, blindly adopting the assumption of homogeneous error terms without exploring possible candidate covariance structures may compromise the ability of an experiment to detect sufficient variation in the response variable. In addition, our approach enhances the accuracy of inferences by providing partitioned analysis of heterogeneous variances and covariance structures, which sometimes are not identical in the data subsets.
Given the increased complexity of research data in various research fields, the application of a linear mixed-model methodology must be in line with data covariance structures for accurate results to be achieved. One of the approaches that has proved to be a reliable tool for managing big data complexity issues is the partitioning approach [12,20,21], in which the traditional homogeneous error variance structure is assumed. The current study extends the new approach to a three-factor treatment structure in a repeated-measures design in which linear mixed-effects models are applicable. In essence, the approach can be extended to cater for repeated-measures experiments in which any number of betweensubjects and within-subjects factors are involved. In most cases, repeated-measures ex-periments do not assume equal and uncorrelated error vectors, since regularly timed measurements taken on the same subject over time are usually correlated [34].
For the fixed-effects partitions, the linear mixed-effects models for a repeated-measures design are fit by the PROC GLM procedure, and a combined analysis of these can be obtained by syncretising the sum of squares and degrees of freedom from the fit models. However, obtaining a combined analysis using the SAS PROC MIXED procedure is impossible using the sum of squares approach, since the PROC MIXED procedure uses a likelihood-based estimation scheme instead of the least-squares method. A comparable alternative to reduction in the sum of squares for the fit model in PROC MIXED is to consider the amount of information retained by the fit model when compared to the null model.

Conclusions
The proposed approach allows for construction and hypothesis testing in repeatedmeasures data if a heterogeneous error structure is assumed. Although the MLE method used in SAS PROC MIXED does not estimate sums of squares, data scrapping based on targeted factor levels work equally well as an alternative approach to obtain the combined analysis. The proposed approach can be adopted as an essential tool for the comparison of new inventions against existing strategies and equipment. It leads to an exploration of the fixed and random levels of the same factor and the subsequent interaction of levels of factors of interest. We can assess the differences between levels of the same factor and understand variation within the same factor. In addition, modelling allows for the assessment of various covariance structures. Although this paper focuses on a limited scope of covariance structures, we pose an open research problem on the application of the proposed approach to other designs of experiments that incorporate more complex and non-diagonal covariance structures.  Data Availability Statement: Data for this research are available in [19,31].

Acknowledgments:
The authors would like to acknowledge the support from the University of South Africa, and from the Sol Plaatje University for the publication costs provided.