Source of data and Variables
The current study’s data comprises 154 adult patients suffering from congestive heart failure, who reported that they had visited a hospital at least thrice within two years. The patients’ pulse rate, respiratory rate and weight were considered as response variables, while age, sex, time, place of residence, New York Heart Association Class, diagnostic history, left ventricle ejection fraction, valvular heart disease, smoking status, diabetes status, diastolic blood pressure and systolic blood pressure were deemed as independent variables.
Statistical Methods of Data Analysis
Statistical analysis contain both descriptive and inferential statistics, such as: summary statistics, Exploratory Data Analysis and model comparison have been used in this study. Joint random effects with LMM have been modeled to infer the joint effect of multivariate longitudinal outcomes of CHF patients.
Separated linear mixed effect model
The model which used in this study was a linear mixed effect model which contains fixed effects, random effects, serial correlation and measurement errors and are referred as β, bi, \({ \varepsilon }_{1i}\)and\({ \varepsilon }_{2i}\) respectively.
The three endpoints are longitudinally measured as a vector of responses,\({Y}_{i}\left(t\right)\) at each occasion. Let \({Y}_{i}=\left[\begin{array}{c}{Y}_{ki}\left(t\right)\\ {Y}_{ji}\left(t\right)\\ {Y}_{li}\left(t\right)\end{array}\right]\), the response vectors for the patient i with \({Y}_{\text{k}\text{i}, } {Y}_{\text{j}\text{i}}, and {Y}_{\text{l}\text{i} }\), the \({Y}_{\text{k}\text{i}, } {Y}_{\text{j}\text{i}}, and {Y}_{\text{l}\text{i} }\)is vector of the end points for k, j, l (k, j, l = 1, 2, 3) with \({n}_{1i}\) =\({n}_{2i}\)=\({n}_{3i}\) = \({n}_{i}\). In generally linear mixed model is written in the form of:
$${Y}_{i}\left(t\right)={X}_{i}{\left(t\right)}^{T}\beta +{Z}_{i}{\left(t\right)}^{T}{b}_{i}+{\epsilon }_{1i}\left(t\right)+{\epsilon }_{2i}\left(t\right) \left(1\right)$$
For each of response variables it can be rewritten in the form of
$${Y}_{ki}\left(t\right)={{\mu }}_{k}\left(\text{t}\right)+{\alpha }_{ki}+{b}_{ki}\left(t\right)+{\epsilon }_{1ki}\left(t\right)+{\epsilon }_{2ki}\left(t\right) \left(2\right)$$
Where\({Y}_{i}\left(t\right)\): Measurement of response in \({i}^{th}\) patient at time t, \({X}_{i}\left(t\right)\): Vector of fixed covariate for \({i}^{th}\) patient at time t,\({Z}_{i}\left(t\right)\) : Vector of random covariate for \({i}^{th}\) patient at time t, β: Vector of unknown parameters associated with fixed covariate (of dimension k), The\({b}_{i}\) is a Vector of unknown parameters associated with random covariate for \({i}^{th}\) patient and ∼N(0,D). It is assumed to be independent of measurement error and serial correlation.
The \({n}_{i}x1\) vector \({\epsilon }_{1i}\) represents serial correlation, the correlation of subsequent measurements within a patient. It is assumed that \({\epsilon }_{1i}\) has a normal distribution with mean zero and covariance matrix \({H}_{i}\) distributed independently of measurement error. The \({n}_{i}x1\) vector \({\epsilon }_{2i}\) stands for the measurement error part, the components of \({\epsilon }_{2i}\)are assumed uncorrelated and normally distributed with mean zero and common variance\({\sigma }^{2}\).
Joint linear mixed effect model
The joint model is extension of separate (univariate) models which used to study the joint evolution of two or more biomarkers. In this study the joint which investigated was bivariate and multivariate linear mixed effect model. For the pairs (k, j) of response variables the bivariate linear mixed effect model is given as:
$$\left.\genfrac{}{}{0pt}{}{{Y}_{ki}\left(t\right)={{\mu }}_{k}\left(\text{t}\right)+{\alpha }_{ki}+{b}_{ki}\left(t\right)+{\epsilon }_{1ki}\left(t\right)+{\epsilon }_{2ki}\left(t\right)}{{Y}_{ji}\left(t\right)={{\mu }}_{j}\left(\text{t}\right)+{\alpha }_{ji}+{b}_{ji}\left(t\right)+{\epsilon }_{1ji}\left(t\right)+{\epsilon }_{2ji}\left(t\right)}\right\} \left(3\right)$$
,\(\text{b}=\left(\genfrac{}{}{0pt}{}{{\alpha }_{i}}{{b}_{i}}\right)=\left(\genfrac{}{}{0pt}{}{\left(\begin{array}{c}{\alpha }_{ki}\\ {\alpha }_{ji}\end{array}\right)}{\left(\begin{array}{c}{b}_{ki}\\ {b}_{ji}\end{array}\right)}\right)\) \(\tilde {N}( 0, { D}_{4x4})\), where D is the variance–covariance matrix for the 2 intercepts and 2 slopes associated with every outcome (covariance matrix for random effects), \(( {\varepsilon }_{1ki} \left(\text{t}\right),{ \varepsilon }_{1ji} \left(\text{t}\right))\) ~ N (0,\({}_{2X2}\)), Where \({}_{2X2}\) is matrix of serial correlation structures, \(J,K for all j\ne k\), \(\left( {\varepsilon }_{2ki} \left(\text{t}\right),{ \varepsilon }_{2ji} \left(\text{t}\right)\right) \tilde N\left(\begin{array}{cc}0& \left[\begin{array}{cc}{{\sigma }^{2}}_{k}& 0\\ 0& {{\sigma }^{2}}_{j}\end{array}\right]\end{array}\right)\).
The multivariate linear mixed effect model is given as:
\(\left.\begin{array}{c}{Y}_{ki}\left(t\right)={{\mu }}_{k}\left(\text{t}\right)+{\alpha }_{ki}+{b}_{ki}\left(t\right)+{\epsilon }_{1ki}\left(t\right)+{\epsilon }_{2ki}\left(t\right)\\ {Y}_{ji}\left(t\right)={{\mu }}_{j}\left(\text{t}\right)+{\alpha }_{ji}+{b}_{ji}\left(t\right)+{\epsilon }_{1ji}\left(t\right)+{\epsilon }_{2ji}\left(t\right)\\ {Y}_{li}\left(t\right)={{\mu }}_{l}\left(\text{t}\right)+{\alpha }_{li}+{b}_{li}\left(t\right)+{\epsilon }_{1li}\left(t\right)+{\epsilon }_{2li}\left(t\right)\end{array}\right\}\)                (4)
,\(\text{b}=\left(\genfrac{}{}{0pt}{}{{\alpha }_{i}}{{b}_{i}}\right)=\left(\genfrac{}{}{0pt}{}{\left(\genfrac{}{}{0pt}{}{{\alpha }_{ki}}{\begin{array}{c}{\alpha }_{ji}\\ {\alpha }_{li}\end{array}}\right)}{\left(\genfrac{}{}{0pt}{}{{b}_{ki}}{\begin{array}{c}{b}_{ji}\\ {b}_{li}\end{array}}\right)}\right)\), \(\tilde \text{N}( 0, { D}_{6x6})\), where D is the variance–covariance matrix for the 3 intercepts and 3 slopes associated with every outcome (covariance matrix for random effects), \(( {\varepsilon }_{1ki} \left(\text{t}\right),{ \varepsilon }_{1ji} \left(\text{t}\right),{ \varepsilon }_{1li}(\text{t}\left)\right)\) ~ N (0,\({}_{3X3}\)), Where \({}_{3X3}\) is matrix of serial correlation structures, \(J,K, and l for all j\ne k\ne l\),
$$\left( {\varepsilon }_{2ki} \left(\text{t}\right),{ \varepsilon }_{2ji} \left(\text{t}\right), { \varepsilon }_{2li} \left(\text{t}\right) \right) \tilde N\left(\begin{array}{cc}0& \left[\begin{array}{ccc}{{\sigma }^{2}}_{k}& 0& 0\\ 0& {{\sigma }^{2}}_{j}& 0\\ 0& 0& {{\sigma }^{2}}_{l}\end{array}\right]\end{array}\right)$$
.
VarianceCovariance and Correlation Structures
In the linear mixed model setup the covariance structure of the repeated measurements can be viewed as consisting of two components that; the structure induced by the random coefficient model, which is typically nonstationary and the residual components which accounts for deviations of the observations about the individual profiles. Most of the time there has been a tendency to use a very simple structure for component (residual components), often independent deviations with common variance. This has implied that the random coefficient model should capture most of the structure of the variance and covariance of the repeated measurements. In this study the researchers relaxed this assumption (independent deviations with common variance) by including the variancecovariance structure for random effects and serial correlation structure to accounts for deviations of the observations about the individual profiles separately. The variance–covariance and correlation matrixes used in this study were listed in the Table 1.
Table 1
The VarianceCovariance/Correlation Structures
Type

Matrix R

Number of Parameters

Unstructured (UN)

UN=\(\left[\begin{array}{ccc}{{\sigma }_{1}}^{2}& {\sigma }_{12}& \dots {\sigma }_{1i}\\ {\sigma }_{21}& {{\sigma }_{2}}^{2}& \dots {\sigma }_{2i}\\ \begin{array}{c}:\\ {\sigma }_{i1}\end{array}& \begin{array}{c}:\\ {\sigma }_{i2}\end{array}& \begin{array}{c}:\\ \dots \end{array} \begin{array}{c}:\\ {{\sigma }_{i}}^{2}\end{array}\end{array}\right]\)

i(i + 1)/2

Compound Symmetry (CS)

CS=\({\sigma }^{2}\left[\begin{array}{ccc}1& \rho & \dots \rho \\ & 1& \dots ⋮\\ & & \ddots \begin{array}{c}\rho \\ 1\end{array}\end{array}\right]\)

2

Heterogeneous Compound symmetry (CSH)

CSH=\(\left[\begin{array}{ccc}{{\sigma }_{1}}^{2}& {\sigma }_{1}{\sigma }_{2}\rho & \dots {\sigma }_{1}{\sigma }_{i}\rho \\ & {{\sigma }_{2}}^{2}& \dots {\sigma }_{2}{\sigma }_{i}\rho \\ & & \begin{array}{cc}\ddots & ⋮\\ & {{\sigma }_{i}}^{2} \end{array} \end{array}\right]\)

i + 1

First Order Autoregressive (AR(1))

AR(1) =\({\sigma }^{2}\left[\begin{array}{ccc}1& \rho & \dots {\rho }^{i1}\\ & 1& \dots {\rho }^{i2}\\ & & \begin{array}{cc}\ddots & ⋮\\ & 1 \end{array} \end{array}\right]\)

2

Toeplitz (TOEP)

TOEP =\(\left[\begin{array}{ccc}{\sigma }^{2}& {\sigma }_{12}& {\sigma }_{13} \dots {\sigma }_{1i}\\ & {\sigma }^{2}& {\sigma }_{12} \dots {\sigma }_{1i1}\\ & & \begin{array}{ccc}{\sigma }^{2} & \dots & ⋮\\ & \ddots & {\sigma }_{12}\\ & & {\sigma }^{2}\end{array}\end{array}\right]\)

T

Variance Components (VC)

VC=\(\left[\begin{array}{ccc}{{\sigma }_{1}}^{2}& 0& \dots 0\\ 0& {{\sigma }_{2}}^{2}& \dots 0\\ \begin{array}{c}:\\ 0\end{array}& \begin{array}{c}:\\ 0\end{array}& \begin{array}{c}:\\ \dots \end{array} \begin{array}{c}:\\ {{\sigma }_{i}}^{2}\end{array}\end{array}\right]\)

Q

Parameter Estimation
To estimate the parameters pertaining to the model used, both maximum likelihood estimation (MLE) and restricted maximum likelihood estimation were performed. The ML method primarily maximizes the log likelihood with respect to the variance parameters, while treating the fixed effects parameters β as constant. After determining the variance parameter estimates, the fixed effects parameters are identified by ascertaining the values of β, which maximize the log likelihood, while treating the variance parameters as constant. Consequently, the variance parameters are estimated by maximizing the REML log likelihood with regard to the variance. Given the nature of the REML likelihood and its treatment of fixed effects as parameters rather than as constants, the resulting variance parameter estimates can be considered as unbiased. Moreover, it reduces the bias inherent in the ML estimates of covariance parameters [14].
Model Comparisons
In the case of repeated or longitudinal data where selecting the best model means not only to select the best mean structure but also the most optimal variance covariance structure for model. In this study the most commonly known model selection criteria which are Akaike Information Criterion (AIC), the Bayesian Information Criterion (BIC), and Loglikelihood ratio test were used to select the best model.