Marginal methods have been widely used for analyzing longitudinal ordinal data due to their simplicity in model assumptions, robustness in inference results, and easiness in the implementation. However, they are often inapplicable in the presence of measurement errors in the variables. Under the setup of longitudinal studies with ordinal responses and covariates subject to misclassification, Chen et al. (2014) developed marginal methods for misclassification adjustments using the second-order estimating equations and proposed a two-stage estimation approach when the validation subsample is available. Parameter estimation is conducted through the Newton-Raphson algorithm, and the asymptotic distribution of the estimators is established. While the methods of Chen et al. (2014) can successfully correct the misclassification effects, its implementation is not accessible to general users due to the lack of a software package. In this paper, we develop an R package, mgee2, to implement the marginal methods proposed by Chen et al. (2014). To evaluate the performance and illustrate the features of the package, we conduct numerical studies.
Analysis of longitudinal ordinal data is a common research topic in health science and survey sampling. Typically, Liang and Zeger (1986) introduced the generalized estimating equations (GEE) method that gave consistent estimation with mild assumptions of the joint distribution of the repeated measurements. This method has been used widely in analyzing longitudinal binary and categorical data. The validity of the GEE method hinges on the critical condition that data are precisely observed, which is commonly infeasible and violated in practice (Yi 2017). Extensive discussions about covariate error (Carroll et al. 2006) and response with binary misclassification (Neuhaus 1999; Chen et al. 2011; Yi 2017) have been conducted in the literature. For example, Neuhaus (1999) investigated the bias due to errors in the response. Yi (2008) proposed a simulation–extrapolation (SIMEX) method to handle both dropout and covariate measurement error problems in longitudinal studies. Furthermore, in Yi (2017 Ch5), the impact of covariate measurement error on longitudinal data analysis was investigated, and methods of addressing covariate measurement error effects were described.
To accommodate effects induced from error-prone correlated ordinal responses and ordinal covariates simultaneously, (Chen et al. 2014) proposed GEE-based methods for the estimation of both mean and association parameters. The proposed methods are based on formulating unbiased second-order estimating functions and solving the resulting equations using the Newton-Raphson algorithm. The asymptotic distributions for the proposed estimators are established. While the methods of (Chen et al. 2014) correct for error effects due to misclassified variables, the methods cannot be used by the analysts without programming the implementation procedures. To expedite the use of the methods for problems in applications, in this paper, we develop an R package, called mgee2, to implement the methods of (Chen et al. 2014).
Our work offers an R package complement to available R packages for analyzing longitudinal data with misclassified observations. It is relevant to but differs from available R packages about measurement error. For example, the package SAMBA, developed by (Beesley and Mukherjee 2020), provides resources for fitting logistic regression with misclassified binary outcomes. The R package misclassGLM implements inferential procedures for generalized linear models with misclassified covariates proposed by (Dlugosz et al. 2017); (Zhang and Yi 2019) developed the package augSIMEX to implement the method proposed by (Yi et al. 2015) for fitting generalized linear models with mixed continuous and discrete covariates subject to mismeasurement.
When the degree of measurement error is very severe, the observed
surrogate measurements are virtually useless, and hence the
corresponding variables may be alternatively treated as subject to
missingness. Regarding the analysis of longitudinal data with missing
observations, packages kml
and kml3d, developed by
(Genolini et al. 2015), describe the implementation procedures of
The article is organized as follows. Section 2 introduces the notations and estimation procedures proposed by (Chen et al. 2014). Section 3 describes the usage of the package mgee2. Section 4 illustrates the package by simulation studies and a real dataset. We finally conclude the article in Section 6.
We first review the notation and formulations of (Chen et al. 2014). For
Let
Let
As a result, the estimating functions for the mean and association
parameters
If the true measurements of the responses and covariates are available,
((3)) and ((4)) can be used for estimation of
Let
If
Case 1 highlights the estimation of
Define
Let
(Chen et al. 2014) established the asymptotic distribution of
We develop an R package, called
mgee2, to implement the
misclassification adjustment method described in the preceding section.
This package requires support from the external packages
MASS (Venables and Ripley 2002),
Matrix (Bates and Maechler 2019), and
ggplot2 (Wickham 2016). Our
mgee2 package mainly
contains two functions, mgee2k
and mgee2v
, respectively,
implementing cases 1 and 2 described in the previous section.
Specifically, mgee2k
implements the method where the misclassification
parameters are given, and mgee2v
implements the misclassification
method for the case where validation data are available to estimate
misclassification probabilities. We now describe the details of these
two functions.
mgee2k
mgee2k
implements the misclassification adjustment method outlined in
Case 1 of the previous section, where the misclassification parameters
are known. In this case, validation data are not required, and only the
observed data of the outcome and covariates are needed for the
implementation.
The function mgee2k
requires the data set to be grouped by the
individual id, id
. The
misclassification matrices for the response and covariate variables are
recorded by the arguments gamMat
and varphiMat
, respectively, which
need to be specified by the user.
To call mgee2k
, we issue the following command,
mgee2k(formula, id, data, corstr="exchangeable", misvariable,
gamMat, varphiMat, maxit=50, tol=1e-3)
where the meaning of each argument is described as follows:
formula
: a formula object which specifies the relationship between
the response and covariates for the observed data.id
: a character object which records individual id in the data.data
: a dataframe or matrix object for the observed data set.corstr
: a character object. The default value is "exchangeable",
corresponding to the structure where the association between two
paired responses is considered to be a constant. The other option is
"log-linear" which indicates the log-linear association between
two paired responses.misvariable
: a character object which names the error-prone
covariate W.maxit
: an integer which specifies the maximum number of
iterations. The default is 50.tol
: a numeric object which indicates the tolerance threshold. The
default is 1e-3.gamMat
: a matrix object which records the misclassification
parameter varphiMat
: a matrix object which records the misclassification
parameter The function mgee2k
returns a list of components:
beta
: the coefficients in the same order as that specified in the
formula for the response and covariates.alpha
: the coefficients for paired responses global odds ratios.
The number of corstr
. When
corstr="exchangeable"
, only one baseline variance
: the variance-covariance matrix of the estimators of all
parameters.convergence
: a logical variable; TRUE if the model converges.iteration
: the number of iterations for the estimates of the model
parameters to converge.call
: an unevaluated function call which consists of the named
function applied to the given arguments.mgee2v
The function mgee2v
does not require the misclassification parameters
to be known, but requires the availability of validation data.
Similar to mgee2k
, the function mgee2v
needs the data set to be
structured by individual id, delta
should be added in the data set, and we
use a column named valid.sample.ind
for this purpose. The column name
of the error-prone covariate misvariable
. To call mgee2v
, we issue the command,
mgee2v(formula, id, data, corstr="exchangeable", misvariable, valid.sample.ind,
y.mcformula, x.mcformula, maxit=50, tol=1e-3)
where the arguments are described as follows:
formula
: a formula object which specifies the relationship between
the response and covariates for the observed data.id
: a character object which records individual id in the data.data
: a dataframe or matrix object for the observed data set.corstr
: a character object. The default value is "exchangeable",
corresponding to the structure where the association between two
paired responses is considered to be a constant. The other option is
"log-linear" which indicates the log-linear association between
two paired responses.misvariable
: a character object which names the error-prone
covariate W.valid.sample.ind
: a string object which names the indicator
variable delta. When a data point belongs to the validation set,
delta = 1; otherwise 0.y.mcformula
: a string object which indicates the misclassification
formula between true response x.mcformula
: a string object which indicates the misclassification
formula between true error-prone covariate X and the surrogate W.maxit
: an integer which specifies the maximum number of
iterations. The default is 50.tol
: a numeric object which indicates the tolerance threshold. The
default is 1e-3.The function mgee2v
returns a list of components:
beta
: the coefficients in the same order as that specified in the
formula for the response and covariates.alpha
: the coefficients for paired responses global odds ratios.
The number of corstr
. When
corstr="exchangeable"
, only one baseline variance
: the variance-covariance matrix of the estimators of all
parameters.convergence
: a logical variable; TRUE if the model converges.iteration
: the number of iterations for the estimates of the model
parameters to converge.call
: an unevaluated function call which consists of the named
function applied to the given arguments.ordGEE2
In addition to developing the package
mgee2 to implement the
methods of (Chen et al. 2014), which accommodate misclassification effects in
inferential procedures, we also implement the naive method of ignoring
the feature of misclassification and call the resulting function
ordGEE2
. This function can be used together with the preceding
described mgee2k
or mgee2v
to evaluate the impact of not addressing
misclassification effects:
ordGEE2(formula, id, data, corstr = "exchangeable", maxit = 50, tol = 0.001)
In this section, we conduct numerical studies to demonstrate the usage
of our developed R package as well as to show supplementary functions
such as summary and plot functions in this package. We first demonstrate
all of the external functions in mgee2
through an example with a
simulated data set, known as obs1
, provided in our package.
The simulated data set, called "obs1
", includes 8 columns and 3000
rows, with each patient having 3 entries of visits. The format of this
data set is as follows.
> head(obs1)
ID Y X treatment visit S W delta
1 1 2 2 1 1 2 2 1
2 1 0 0 1 2 0 0 1
3 1 <NA> <NA> 1 3 1 2 0
4 2 <NA> <NA> 1 1 1 0 0
5 2 <NA> <NA> 1 2 0 1 0
6 2 <NA> <NA> 1 3 0 0 0
> summary(obs1)
ID Y X treatment visit
Min. : 1.0 0 : 352 0 : 444 0:1500 1:1000
1st Qu.: 250.8 1 : 283 1 : 269 1:1500 2:1000
Median : 500.5 2 : 256 2 : 178 3:1000
Mean : 500.5 NA's:2109 NA's:2109
3rd Qu.: 750.2
Max. :1000.0
S W delta
0:1181 0:1460 Min. :0.000
1: 955 1: 944 1st Qu.:0.000
2: 864 2: 596 Median :0.000
Mean :0.297
3rd Qu.:1.000
Max. :1.000
Here, delta
is 1 when the subject is in the validation set and 0 otherwise.
About
In the case corstr = "exchangeable"
, the association, defined as in
((2)), between paired responses is assumed to be
corstr = "log-linear"
, the association is assumed to
be
We now apply mgee2k
and mgee2v
, in contrast to ordGEE2
, to fit the
data to the models, respectively. The results are displayed as follows.
In the summary tables for the R output, we use "Y>=1
" and "Y>=2
"
to denote the coefficients Delta
" correspond to the parameter
To use function mgee2k
, we need to specify the misclassification
matrices beforehand. Here, we set the misclassification matrices the
same as used in the simulation process.
> data(obs1)
> obs1$visit <- as.factor(obs1$visit)
> obs1$treatment <- as.factor(obs1$treatment)
> obs1$S <- as.factor(obs1$S)
> obs1$W <- as.factor(obs1$W)
> ## set misclassification parameters to be known.
> varphiMat <- gamMat <- log( cbind(0.04/0.95, 0.01/0.95,
+ 0.95/0.03, 0.02/0.03,
+ 0.04/0.01, 0.95/0.01) )
> mgee2k.fit = mgee2k(formula = S~W+treatment+visit, id = "ID", data = obs1,
+ corstr = "exchangeable", misvariable = "W", gamMat = gamMat,
+ varphiMat = varphiMat)
> summary(mgee2k.fit)
Call:
mgee2k(formula = S ~ W + treatment + visit, id = "ID", data = obs1,
corstr = "exchangeable", misvariable = "W", gamMat = gamMat,
varphiMat = varphiMat)
Summary table of the estimation
Estimate Std.Err Z value Pr(>z)
Y>=1 0.70889 0.08591 8.251 2.22e-16 ***
Y>=2 -0.67521 0.08625 -7.828 4.88e-15 ***
W1 0.58667 0.08719 6.729 1.71e-11 ***
W2 0.94948 0.09745 9.743 < 2e-16 ***
treatment1 -0.70554 0.09114 -7.742 9.77e-15 ***
visit2 -0.24147 0.07735 -3.122 0.0018 **
visit3 -0.62480 0.07571 -8.253 2.22e-16 ***
Delta 1.22606 0.12231 10.024 < 2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
To use mgee2v
, a column of indicator variable should be specified in
valid.sample.ind
.
> data(obs1)
> obs1$visit <- as.factor(obs1$visit)
> obs1$treatment <- as.factor(obs1$treatment)
> obs1$S <- as.factor(obs1$S)
> obs1$W <- as.factor(obs1$W)
> mgee2v.fit = mgee2v(formula = S~W+treatment+visit, id = "ID", data = obs1,
+ y.mcformula = "S~1", x.mcformula = "W~1",
+ misvariable = "W", valid.sample.ind = "delta",
+ corstr = "exchangeable")
> summary(mgee2v.fit)
Call:
mgee2v(formula = S ~ W + treatment + visit, id = "ID",
data = obs1, corstr = "exchangeable", misvariable = "W",
valid.sample.ind = "delta", y.mcformula = "S~1", x.mcformula = "W~1")
Summary table of the estimation
Estimate Std.Err Z value Pr(>z)
Y>=1 0.64876 0.08851 7.330 2.30e-13 ***
Y>=2 -0.68226 0.08703 -7.839 4.44e-15 ***
W1 0.56507 0.08140 6.942 3.88e-12 ***
W2 0.98411 0.09305 10.577 < 2e-16 ***
treatment1 -0.68153 0.09052 -7.529 5.11e-14 ***
visit2 -0.24694 0.07483 -3.300 0.000966 ***
visit3 -0.60027 0.07335 -8.184 2.22e-16 ***
Delta 1.22862 0.12160 10.103 < 2e-16 ***
---
Signif. codes:
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> naigee.fit = ordGEE2(formula = S~W+treatment+visit, id = "ID",
+ data = obs1, corstr = "exchangeable")
> summary(naigee.fit)
Call:
ordGEE2(formula = S ~ W + treatment + visit, id = "ID",
data = obs1, corstr = "exchangeable")
Summary table of the estimation
Estimate Std.Err Z value Pr(>z)
Y>=1 0.73276 0.07990 9.171 < 2e-16 ***
Y>=2 -0.69330 0.08004 -8.662 < 2e-16 ***
W1 0.51237 0.07354 6.967 3.23e-12 ***
W2 0.84890 0.08582 9.892 < 2e-16 ***
treatment1 -0.65954 0.08511 -7.749 9.33e-15 ***
visit2 -0.22766 0.07241 -3.144 0.00167 **
visit3 -0.58407 0.07052 -8.282 2.22e-16 ***
Delta 1.06616 0.09846 10.828 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 '' 1
plot_model
We use the function plot_model
to compare the results obtained from
the three functions:
> plot_model(naigee.fit)
> plot_model(mgee2.fit)
> plot_model(mgee2v.fit)
It is helpful to compare the odds ratios when there are multiple
covariates. We use the function plot_model
to visualize the odds
ratios. The estimated odds ratios for this simulated data set across the
three methods are displayed in Figure 1, 2,
and 3. The red dot gives the odds ratio of each
covariate. The horizontal blue line measures the length of each
confidence interval. The vertical axes of the graphs indicate the
descending order of the covariates. In other words, the red points from
the lowest to the highest in the graph represent the first covariate,
the second covariate, and so on. It is seen that the three methods yield
similar odds ratios.
To further compare the three methods, a simulation study is conducted.
We run 500 simulations where each data set includes 1000 subjects, with
three visits for each subjects. obs1
is one example of the simulated
data. The true values of the coefficients are reported in
1:
2 reports the simulation results for the care with a 5%
misclassification rate set for both the response and covariate
variables, where Bias% records a bias in percentage, EV represents an
empirical variance, AMV stands for an average of model-based variance,
and CR records a coverage rate of 95% confidence intervals. Simulation
results show that the mgee2
and mgee2v
perform better than the naive
method ordGEE2
, and they produce reasonable results.
ordGEE2 | mgee2 | mgee2v | ||||||||||
Bias% | EV | AMV | CR | Bias% | EV | AMV | CR | Bias% | EV | AMV | CR | |
3.119 | 0.007 | 0.007 | 0.942 | -0.915 | 0.008 | 0.008 | 0.944 | -2.223 | 0.008 | 0.008 | 0.951 | |
3.226 | 0.007 | 0.007 | 0.946 | 1.562 | 0.008 | 0.008 | 0.940 | 3.556 | 0.009 | 0.008 | 0.947 | |
-12.112 | 0.006 | 0.006 | 0.784 | 1.238 | 0.008 | 0.008 | 0.942 | -6.933 | 0.024 | 0.014 | 0.924 | |
-9.810 | 0.007 | 0.008 | 0.754 | 1.433 | 0.009 | 0.010 | 0.964 | 3.393 | 0.016 | 0.011 | 0.941 | |
-7.032 | 0.008 | 0.007 | 0.922 | -0.456 | 0.009 | 0.008 | 0.954 | -0.138 | 0.009 | 0.008 | 0.949 | |
-6.311 | 0.006 | 0.005 | 0.932 | 0.071 | 0.006 | 0.006 | 0.938 | -0.364 | 0.006 | 0.006 | 0.932 | |
-6.630 | 0.005 | 0.005 | 0.908 | 0.056 | 0.006 | 0.006 | 0.964 | 0.143 | 0.006 | 0.006 | 0.962 | |
-13.130 | 0.009 | 0.010 | 0.690 | 0.217 | 0.014 | 0.015 | 0.954 | 1.257 | 0.016 | 0.017 | 0.956 |
In addition to the preceding simulation with a misclassification rate of
5%, we conducted another simulation with the same parameters except that
the misclassification rate is changed to be 20%, and
corstr = "log-linear"
. The results are reported in Table
3, which shows more noticeable differences in
implementing the three functions, ‘ordGEE2’, ‘mgee2’, and ‘mgee2v’.
ordGEE2 | mgee2 | mgee2v | ||||||||||
Bias% | EV | AMV | CR | Bias% | EV | AMV | CR | Bias% | EV | AMV | CR | |
9.589 | 0.007 | 0.007 | 0.866 | 0.748 | 0.015 | 0.015 | 0.952 | 0.872 | 0.015 | 0.016 | 0.966 | |
11.131 | 0.006 | 0.007 | 0.842 | -0.210 | 0.014 | 0.015 | 0.958 | -0.523 | 0.015 | 0.016 | 0.958 | |
-48.891 | 0.005 | 0.006 | 0.000 | -1.233 | 0.027 | 0.029 | 0.964 | -0.874 | 0.023 | 0.023 | 0.940 | |
-43.506 | 0.008 | 0.008 | 0.000 | -0.400 | 0.026 | 0.027 | 0.958 | -0.399 | 0.023 | 0.023 | 0.946 | |
-26.284 | 0.006 | 0.006 | 0.346 | 0.364 | 0.011 | 0.012 | 0.960 | 0.084 | 0.011 | 0.011 | 0.940 | |
-24.870 | 0.006 | 0.006 | 0.832 | 1.703 | 0.012 | 0.011 | 0.938 | 1.396 | 0.010 | 0.010 | 0.948 | |
-26.893 | 0.006 | 0.006 | 0.314 | 0.228 | 0.011 | 0.012 | 0.954 | 0.144 | 0.010 | 0.011 | 0.968 | |
-53.210 | 0.009 | 0.009 | 0.000 | 1.873 | 0.078 | 0.068 | 0.942 | 1.034 | 0.052 | 0.066 | 0.976 | |
-73.139 | 0.006 | 0.006 | 0.042 | -0.353 | 0.052 | 0.047 | 0.948 | -1.241 | 0.037 | 0.049 | 0.970 | |
-59.955 | 0.011 | 0.011 | 0.000 | 0.256 | 0.075 | 0.075 | 0.944 | 0.188 | 0.053 | 0.079 | 0.978 |
To illustrate the usage of the developed R package, we analyze a dataset arising from the Framingham Heart Study, obtained from the NIH website (https://biolincc.nhlbi.nih.gov/teaching/). Similar to (Chen et al. 2014), we consider those 915 male patients who completed both exams #2 and #3, and age between 31 and 65 at the entry of the study. The response variable, HBP, is a categorical variable indicating the status of the systolic blood pressure (SBP), where HBP=0 if SBP is below 140 mmHg, HBP=1 if SBP is between 140 mmHg and 159 mmHg, and HBP=2 if SBP is larger than 160 mmHg.
We are interested in understanding the relationship between HBP and covariates, including the serum cholesterol level (CHOL), age, and the current smoking status (CURSMOKE), as well as the examination status, denoted as "Exam3". CHOL is classified as three categories, with 0, 1, and 2 representing normal (less than 200 mg/dL), borderline high (200-239mg/dL), and hypercholesterolemia (greater than 240 mg/dL), respectively. Exam3 is a dummy variable, with 1 indicating observations for exam 3 and 0 for exam 2.
First, we visualize how SBP may change with age by stratifying the study subjects into different categories according to the exam time, smoking status, or CHOL. To see the trend, we display simple linear regression lines that fit scattered points of SBP against AGE for patients in each category, as shown in Figure 4. Except for the patients with CHOL value 0 and CURSMOKE value 0 at exam 2, there is generally an upward trend of SBP versus age for each category, though the degree varies. While each patient takes 2 exams, the time interval between two exams is different from patient to patient. To reflect this feature of different gap times for the study subjects, in Figure 5 we further display spaghetti plots (Hedeker and Gibbons 2006) for patients in different categories, where the two endpoints of each black line segment mark SBP and age for the corresponding study subject at exams 2 and 3 in each category, respectively. The blue curve represents the loess smooth curve in each panel to show the trend of SBP against AGE. The loess smooth function is a tool to create smooth lines for scattered plots using polynomial approximations. The code for producing Figures 4 and 5 is included in the help file of data set heart in our R package.
Next, we use the proportional odds model to examine how SBP may be
quantitatively associated with the covariates. For the
As defined in ((1)), consider the model
ordGEE2 | mgee2k | ||||||||
Est. | SD | p-vlue | Est. | SD | p-vlue | ||||
-4.291 | 0.635 | <0.001 | -4.943 | 0.737 | <0.001 | ||||
-5.623 | 0.638 | <0.001 | -6.195 | 0.740 | <0.001 | ||||
0.068 | 0.133 | 0.608 | 0.117 | 0.180 | 0.515 | ||||
0.352 | 0.140 | 0.012 | 0.474 | 0.186 | 0.011 | ||||
0.063 | 0.012 | <0.001 | 0.071 | 0.014 | <0.001 | ||||
-0.044 | 0.105 | 0.673 | -0.042 | 0.118 | 0.722 | ||||
0.145 | 0.097 | 0.133 | 0.182 | 0.110 | 0.097 | ||||
2.301 | 0.207 | <0.001 | 2.301 | 0.318 | <0.001 |
The data set used in our example is included in our package called
"heart
". To demonstrate the usage of the developed package, we
perceive that the response HBP level and the covariate cholesterol level
are prone to misclassification. Since this example does not have a
validation data set, we only analyze the data using the naive method,
"ordGEE2
", and the corrected method with a specified known
misclassification rate, "mgee2k
", where the misclassification rates
for both the outcome and the covariate are assumed to be 5%, and the
exchangeable dependence structure is considered. The analysis results
are shown in Table 4. Overall, the naive method and the
corrected method indicate the same significant health factors, yet the
magnitude of the coefficient estimates and their standard errors are
different. Higher cholesterol levels and older ages appear to be
positively correlated with high blood pressure.
Analysis of longitudinal ordinal data is important for research in health science, epidemiological studies, and social science. Marginal analysis using generalized estimating equations has been extensively employed in applications. However, such a strategy is challenged by the presence of mismeasurement of variables. To address this challenge, (Chen et al. 2014) developed estimation methods for analyzing correlated ordinal responses and ordinal covariates, which are subject to misclassification.
To allow analysts to apply the useful methods of (Chen et al. 2014) without doing individual codes, we develop an R package mgee2 to implement the methods for general use. Our package provides three methods for estimation, including the two methods of corrections for misclassification effects, as opposed to the naive method, which disregards the feature of mismeasurement in variables. The package can be used for modeling longitudinal ordinal data with misclassified response and covariates. It provides consistent estimation results by directly inputting the data under required assumptions.
The authors thank the review team for the helpful comments on the initial version. Yi’s research was supported by the Natural Sciences and Engineering Research Council of Canada (NSERC). Yi is Canada Research Chair in Data Science (Tier 1). Her research was undertaken, in part, thanks to funding from the Canada Research Chairs program.
mgee2, SAMBA, misclassGLM, augSIMEX, kml, kml3d, gee, wgeesel, swgee, MASS, Matrix, ggplot2, mgee2k, mgee2v, ordGEE2
Cluster, Distributions, Econometrics, Environmetrics, MissingData, MixedModels, NumericalMathematics, Phylogenetics, Psychometrics, Robust, Spatial, TeachingStatistics
This article is converted from a Legacy LaTeX article using the texor package. The pdf version is the official version. To report a problem with the html, refer to CONTRIBUTE on the R Journal homepage.
Text and figures are licensed under Creative Commons Attribution CC BY 4.0. The figures that have been reused from other sources don't fall under this license and can be recognized by a note in their caption: "Figure from ...".
For attribution, please cite this work as
Xu, et al., "mgee2: An R package for marginal analysis of longitudinal ordinal data with misclassified responses and covariates", The R Journal, 2021
BibTeX citation
@article{RJ-2021-093, author = {Xu, Yuliang and Liu, Shuo Shuo and Yi, Grace Y.}, title = {mgee2: An R package for marginal analysis of longitudinal ordinal data with misclassified responses and covariates}, journal = {The R Journal}, year = {2021}, note = {https://rjournal.github.io/}, volume = {13}, issue = {2}, issn = {2073-4859}, pages = {471-484} }