Access the full text.
Sign up today, get DeepDyve free for 14 days.
E. Ziegel (2002)
Generalized Linear ModelsTechnometrics, 44
Atanu Bhattacharjee (2020)
Longitudinal Data AnalysisBayesian Approaches in Oncology Using R and OpenBUGS
Xin Wei (2012)
%PROC_R: A SAS Macro that Enables Native R Programming in the Base SAS EnvironmentJournal of Statistical Software, 046
A. Ziegler, G. Arminger (1995)
Analyzing the Employment Status with Panel Data from the GSOEP: A Comparison of the MECOSA and the GEE1 Approach for Marginal Models, 64
F ACITO (1986)
111Journal of Marketing Research, 23
L. Wei, D. Stram (1988)
Analysing repeated measurements with possibly missing observations by modelling marginal distributions.Statistics in medicine, 7 1-2
(2009)
Ideal Point Discriminant Analysis with a Special Emphasis on Visualization
A SOMMER, J KATZ, I TARWOTJO (1984)
Increased Risk of Respiratory Disease and Diarrhea in Children with Preexsting Mild Vitamin A DeficiencyAmerican Society for Clinical Nutrition, 40
Timo Oertzen, C. Hertzog, U. Lindenberger, P. Ghisletta (2010)
The effect of multiple indicators on the power to detect inter-individual differences in change.The British journal of mathematical and statistical psychology, 63 Pt 3
S. Bull (1998)
Regression models for multiple outcomes in large epidemiologic studies.Statistics in medicine, 17 19
A. Gifi (1990)
NONLINEAR MULTIVARIATE ANALYSIS
A. Boomsma, J. Hoogland (2001)
The robustness of LISREL modeling revisted.
(1992)
Revised NEO Personality Inventory (NEO-PRI) and NEO Five-Factor Inventory (NEO- FFI) Professional Manual, Odessa, FL: Psychological Assessment Resources
I. Plewis (1996)
Statistical methods for understanding cognitive growth: A review, a synthesis and an applicationBritish Journal of Mathematical and Statistical Psychology, 49
C. Braak, P. Verdonschot (1995)
Canonical correspondence analysis and related multivariate methods in aquatic ecologyAquatic Sciences, 57
M. Rooij, H. Worku (2012)
A warning concerning the estimation of multinomial logistic models with correlated responses in SASComputer methods and programs in biomedicine, 107 2
Søren Højsgaard, U. Halekoh, Jun Yan (2005)
The R Package geepack for Generalized Estimating EquationsJournal of Statistical Software, 15
F. Acito, Ronald Anderson (1986)
A Simulation Study of Factor Score IndeterminacyJournal of Marketing Research, 23
M. Sherman, S. Cessie (1997)
A comparison between bootstrap methods and generalized estimating equations for correlated outcomes in generalized linear modelsCommunications in Statistics - Simulation and Computation, 26
Guang Cheng, Zhuqing Yu, Jianhua Huang (2013)
The cluster bootstrap consistency in generalized estimating equationsJ. Multivar. Anal., 115
P. Spinhoven, E. Penelo, M. Rooij, B. Penninx, J. Ormel (2013)
Reciprocal effects of stable and temporary components of neuroticism and affective disorders: results of a longitudinal cohort studyPsychological Medicine, 44
W. Pan (2001)
Akaike's Information Criterion in Generalized Estimating EquationsBiometrics, 57
D. Elliott, D. Huizinga, S. Menard (1991)
Multiple Problem Youth: Delinquency, Substance Use, and Mental Health Problems
KY LIANG, SL ZEGER (1986)
Longitudinal Data Analysis Using Generalised Linear ModelsBiometrika, 73
(1978)
Multidimensional Scaling, Sage Publications
A AGRESTI (2002)
Categorical Data Analysis
K. Liang, S. Zeger, B. Qaqish (1992)
Multivariate Regression Analyses for Categorical DataJournal of the royal statistical society series b-methodological, 54
T. Park (1994)
Multivariate regression models for discrete and continuous repeated measurementsCommunications in Statistics-theory and Methods, 23
Alfred Sommer, J. Katz, I. Tarwotjo (1984)
Increased risk of respiratory disease and diarrhea in children with preexisting mild vitamin A deficiency.The American journal of clinical nutrition, 40 5
PGM HEIJDEN, A MOOIJAART, Y TAKANE (1994)
Correspondence Analysis in the Social Sciences
(1971)
The Biplot Graphical Display of Matrices with Application to Principal Component Analysis
H AKAIKE (1973)
Proceedings of the Second International Symposium on Information Theory
H. Chaiklin (1992)
Multiple Problem Youth: Delinquency, Substance Use, and Mental Health ProblemsJournal of Nervous and Mental Disease, 180
A BOOMSMA, JJ HOOGLAND (2001)
Structural Equation Modeling: Present and Future
R. Krueger (1999)
The structure of common mental disorders.Archives of general psychiatry, 56 10
Özgür Asar, Özlem Ilk (2013)
mmm: An R package for analyzing multivariate longitudinal data with multivariate marginal modelsComputer methods and programs in biomedicine, 112 3
(2009)
The Structure of CommonMental Disorders: A Replication Study in a Community Sample of Adolescents and Young Adults
A. Ziegler, C. Kastner, M. Blettner (1998)
The Generalised Estimating Equations: An Annotated BibliographyBiometrical Journal, 40
S. Lipsitz, Kyungeun Kim, L. Zhao (1994)
Analysis of repeated categorical data using generalized estimating equations.Statistics in medicine, 13 11
P. Spinhoven, M. Rooij, W. Heiser, J. Smit, B. Penninx (2009)
The role of personality in comorbidity among anxiety and depressive disorders in primary care and specialty care: a cross-sectional analysis.General hospital psychiatry, 31 5
PT COSTA, RR MCCRAE (1992)
Revised NEO Personality Inventory (NEO-PRI) and NEO Five-Factor Inventory (NEO- FFI) Professional Manual
A. Hubbard, J. Ahern, N. Fleischer, M. Laan, Sheri Lippman, N. Jewell, T. Bruckner, W. Satariano (2010)
To GEE or Not to GEE: Comparing Population Average and Mixed Models for Estimating the Associations Between Neighborhood Risk Factors and HealthEpidemiology, 21
E. Chao (2003)
Generalized Estimating EquationsTechnometrics, 45
H. Akaike (1973)
Information Theory and an Extension of the Maximum Likelihood Principle, 1
R. Team (2014)
R: A language and environment for statistical computing.MSOR connections, 1
(2011)
Understanding Biplots
B. Penninx, A. Beekman, J. Smit, F. Zitman, W. Nolen, P. Spinhoven, P. Cuijpers, P. Jong, H. Marwijk, W. Assendelft, K. Meer, P. Verhaak, M. Wensing, R. Graaf, W. Hoogendijk, J. Ormel, R. Dyck (2008)
The Netherlands Study of Depression and Anxiety (NESDA): rationale, objectives and methodsInternational Journal of Methods in Psychiatric Research, 17
C. Braak (1986)
Canonical Correspondence Analysis: A New Eigenvector Technique for Multivariate Direct Gradient AnalysisEcology, 67
H. Worku, M. Rooij (2017)
Properties of Ideal Point Classification Models for Bivariate Binary DataPsychometrika, 82
I BORG, PJF GROENEN (2005)
Modern Multidimensional Scaling: Theory and Applications
JC GOWER, DJ HAND (1996)
Biplots
(1994)
Correspondence Analysis and Contingency Models
Journal of Classiﬁcation 35:124-146 (2018) DOI: 10.1007/s00357-018-9251-4 A Multivariate Logistic Distance Model for the Analysis of Multiple Binary Responses Hailemichael M. Worku Leiden University, The Netherlands Mark de Rooij Leiden University, The Netherlands Abstract: We propose a Multivariate Logistic Distance (MLD) model for the analy- sis of multiple binary responses in the presence of predictors. The MLD model can be used to simultaneously assess the dimensional/factorial structure of the data and to study the effect of the predictor variables on each of the response variables. To en- hance interpretation, the results of the proposed model can be graphically represented in a biplot, showing predictor variable axes, the categories of the response variables and the subjects’ positions. The interpretation of the biplot uses a distance rule. The MLD model belongs to the family of marginal models for multivariate responses, as opposed to latent variable models and conditionally speciﬁed models. By setting the distance between the two categories of every response variable to be equal, the MLD model becomes equivalent to a marginal model for multivariate binary data estimated using a GEE method. In that case the MLD model can be ﬁtted using existing sta- tistical packages with a GEE procedure, e.g., the genmod procedure from SAS or the geepack package from R. Without the equality constraint, the MLD model is a general model which can be ﬁtted by its own right. We applied the proposed model to empirical data to illustrate its advantages. Keywords: Multivariate binary data; Biplots; Multidimensional scaling; Multidi- mensional unfolding; Marginal model; Clustered bootstrap; Generalized estimating equations. Corresponding Author’s Address: H.M. Worku, Psychology Institute, Methodology and Statistics Unit, Leiden University, Wassenaarseweg 52, Box 9555, 2300 RB, Leiden, The Netherlands, e-mail: hmetiku@yahoo.com. P ublished online: 26 March 2018 125 Multivariate Logistic Distance Models 1. Introduction Multivariate binary data with multiple binary response variables and one or more predictor variables, are often collected in empirical sciences such as psychology, criminology, epidemiology, life sciences and medicine. In the Netherlands Study of Depression and Anxiety (NESDA), for example, data were collected to investigate the interplay between personality traits and co-morbidity of depressive and anxiety disorders (Pennix et al., 2008; Spin- hoven et al., 2009). Another study in which multivariate binary data arises is the Indonesian Children’s Study (ICS: Sommer, Katz, and Tarwotjo, 1984; Liang, Zeger, and Qaqish, 1992) where over three-thousand children were medically examined to investigate whether they had respiratory infection, diarrhoeal infection, and xerophthalmia. The aim of the ICS study was to investigate whether vitamin A deﬁciency places children at increased risk of respiratory and diarrhoeal infections. The availability of the multivariate normal distribution for multivari- ate interval responses, makes application of maximum likelihood-based sta- tistical models on such data relatively easy. However, for binary responses, no multivariate distribution is available and therefore estimation becomes more difﬁcult. Liang and Zeger (1986) proposed Generalized Estimating Equations (GEE) for marginal modelling of correlated categorical data. GEE is a quasi-likelihood (QL) estimation method that does not require speciﬁca- tion of a particular multivariate distribution. It is widely used as a standard approach for ﬁtting marginal models on multivariate data (Ziegler, Kastner, and Blettner, 1998; Fitzmaurice et al., 2008; Ziegler, 2011). The GEE ap- proach, however, does not allow for a dimensional approach to analysis. Often researchers have theories how different response variables belong to one underlying dimension, factor, or latent variable. For the dimensional approach often latent variable models are used, such as structural equation models or item response models. These mod- els explicitly deﬁne underlying dimensions. However, these models make distributional assumptions of the latent dimensions or assume an underlying distribution for the dichotomous responses or both. Such assumptions are often unveriﬁable, i.e. we cannot check the assumptions using the data. In this paper, we will develop a dimensional model for multivariate bi- nary data within the marginal framework. The model does not make unver- iﬁable assumptions. The model will be developed within a distance frame- work, but we show it can also be seen as a speciﬁc marginal model. To enhance interpretation, a biplot is developed to accompany the model that visualizes the result. De Rooij (2009) proposed the Ideal Point Classiﬁcation (IPC) model for analyzing a multinomial response variable in the presence of predic- 126 H.M. Worku and M. de Rooij tors. The IPC is a probabilistic distance model based on a two-mode dis- tance function. De Rooij (2009) also showed that a simple logistic regres- sion for binary response variable can be written as a unidimensional IPC model. Worku and De Rooij (2016) extended the IPC model to the analysis of two binary response variables, i.e., the bivariate, binary data setting, and showed that a new parameterization of the IPC model recovered both the marginal probabilities and the association structure of bivariate binary data well. However, this parameterization cannot be easily extended to handling multivariate binary data because all the possible pairwise and higher order association terms must be speciﬁed in the likelihood function, which makes the model complex and therefore hard to estimate. Therefore, in this paper we propose a Multivariate Logistic Distance (MLD) model for analyzing multivariate binary data that extends marginal models for multivariate data. The MLD model uniﬁes two domains of sta- tistical methods, i.e., Multidimensional Scaling (MDS: Kruskal and Wish, 1978; Borg and Groenen, 2005) and Generalized Linear Model (GLM: McCullagh and Nelder, 1989; Agresti, 2002). As a form of regularization, the MLD model allows for dimension reduction and therefore less parame- ters are estimated compared to the existing marginal models for multivari- ate data. Moreover, the model enhances interpretation by using a biplot (Gabriel, 1971; Gower and Hand, 1996; Gower, Lubbe, and Le Roux, 2011) based on a distance interpretation. Unlike existing marginal models for multivariate data, the MLD model can be used for assessing the factorial/dimensional structure of multivari- ate data. In the area of mental disorders (with the NESDA data as exam- ple), clinical psychologists and epidemiologists are often interested in co- morbidity and how co-morbidity is related to risk factors such as personality traits (Krueger, 1999; Beesdo-Baum et al., 2009; Spinhoven et al., 2013). Three candidate theories about the co-morbidity of mental disorders have been proposed, i.e., (1) a 2-dimensional structure with one dimension rep- resenting distress and the other one fear (d/f); (2) a different 2-dimensional structure with one dimension representing depression and the other one anx- iety (d/a); and (3) an unidimensional structure where all the disorders are represented by a single dimension. The MLD model can be used to repre- sent such theories within a uniﬁed framework, i.e., the candidate theories can be compared using appropriate statistics, and at the same time the MLD model allows for a direct relationship between co-morbidity of mental dis- orders and the predictor variables. The paper is organized as follows. Section 2 develops the multivari- ate logistic distance model, investigates the link with marginal model for multivariate binary data estimated using a GEE method, and discusses the construction of biplots for the multivariate logistic model. In Section 3, the 127 Multivariate Logistic Distance Models proposed model is ﬁtted to empirical data and the results are interpreted us- ing the estimated parameters and a graphical representation. We conclude in Section 4 with a discussion. 2. Multivariate Logistic Regression in a Distance Framework 2.1 Logistic Regression as a Distance Model Logistic regression is a standard method for modelling dichotomous response data. Let y denote the observed value for a binary dependent vari- able Y for subject i,where i =1, 2,... ,N. Logistic regression models the probability of a category conditional on the value of a predictor variable x , Pr(y =1|x )= π(x ),i.e., i i i ∗ ∗ exp(β + β x ) 0 1 π(x )= , (1) ∗ ∗ 1+ exp(β + β x ) 0 1 ∗ ∗ where β and β are the intercept and the regression coefﬁcient, respectively. 0 1 Logistic regression can easily be generalized to accommodate multiple pre- T ∗ T dictors, x =(x ,x ,... ,x ) , and thus π(x )= exp(β + x β )/(1 + i i1 i2 ip i ∗ ∗ ∗ T exp(β + x β )),where β is now a vector with regression coefﬁcients. De Rooij (2009) showed that logistic regression can be expressed as a distance model in a joint space with points representing the two categories of the response variable and points representing the subjects. In this sec- tion, we revisit this relationship and in Section 2.2 discuss an extension for multivariate binary responses. Let us deﬁne a joint unidimensional space for subjects and the classes of the response variables. Denote by η the coordinate of the position for subject i and by γ the coordinate of the position of one category and by γ 0 1 the coordinate of the position of the other category of the binary response variable. Deﬁne δ and δ to be the squared Euclidean distances between i0 i1 the position of subject i and the two categories respectively. That is, δ =(η − γ ) ; i1 i 1 (2) δ =(η − γ ) . i0 i 0 With these two distances we can deﬁne the following probability model exp(−0.5δ ) i1 π(x )= . (3) exp(−0.5δ )+exp(−0.5δ ) i0 i1 The smaller the relative distance between a person point and a class point, the larger the probability that the subject belongs to that class. Therefore, 128 H.M. Worku and M. de Rooij the class probability is inversely related to the squared Euclidean distance between the points. The coordinate for subject i, η , is assumed to be a linear combination of the predictor variable x ,i.e., η = β + β x or in case of multiple i i 0 1 i predictors η = β + x β. The parameters of the distance model are the i 0 regression weights and the category points. An important tool in the interpretation of probability models is the log-odds. The log-odds representation of the distance model becomes, π(x ) log =0.5δ − 0.5δ i0 i1 1 − π(x ) 2 2 = η (γ − γ )+ 0.5(γ − γ ) i 1 0 (4) 0 1 2 2 =(β + β x )(γ − γ )+0.5(γ − γ ) 0 1 i 1 0 0 1 2 2 = β (γ − γ )+0.5(γ − γ )+ β (γ − γ )x . 0 1 0 1 1 0 i 0 1 In the case of multiple predictors, the logistic distance model takes the same form, having an intercept and extra slopes for the additional predictors. For example, with two predictors x =(x ,x ) , the distance model becomes, i i1 i2 π(x ) 2 2 log = β (γ − γ )+0.5(γ − γ ) 0 1 0 0 1 1 − π(x ) (5) + β (γ − γ )x + β (γ − γ )x . 1 1 0 i1 2 1 0 i2 For a unit increase in x , the log-odds in the distance model changes by i1 ∗ 2 2 β (γ − γ ), similarly for x . By setting β = β (γ − γ )+0.5(γ − γ ) 1 1 0 i2 0 1 0 0 0 1 and β = β (γ − γ ), a standard logistic regression is obtained. 1 1 0 The logistic distance model (4) is not identiﬁed and therefore an iden- tiﬁability constraint must be imposed. For example, with β =2 and ∗ ∗ (γ − γ )=1, β =2.The same value β =2 can also be obtained 1 0 1 1 when β =0.5 and (γ − γ )= 2. By imposing an identiﬁability con- 1 1 0 straint on the class points, the logistic distance model can be identiﬁed, for example by setting γ =1 and γ =0. The logistic distance model is now 1 1 identiﬁed and its relationship with the univariate logistic model presented in (1) becomes β = β − 0.5; (6) β = β . 2.2 Multivariate Extension of the Distance Model In this section, the logistic distance model for a single response vari- able will be extended to handling multivariate binary data. Suppose y = (y ,y ,... ,y ,... ,y ) denotes the multivariate responses observed on i1 i2 ij iJ 129 Multivariate Logistic Distance Models Table 1. The structure of multivariate data in long format. Predictor variables SID Index Response x x x 1 2 p 1 R y x x ... x 1 11 11 12 1p 1 R y x x ... x 2 12 11 12 1p 1 R y x x ... x 3 13 11 12 1p 1 R y x x ... x 4 14 11 12 1p 1 R y x x ... x 5 15 11 12 1p . . . . . . . . . . . . . . . . . . . . . i R y x x ... x 1 i1 i1 i2 ip i R y x x ... x 2 i2 i1 i2 ip i R y x x ... x 3 i3 i1 i2 ip i R y x x ... x 4 i4 i1 i2 ip i R y x x ... x 5 i5 i1 i2 ip . . . . . . . . . . . . . . . . . . . . . n R y x x ... x 1 n1 n1 n2 np n R y x x ... x 2 n2 n1 n2 np n R y x x ... x 3 n3 n1 n2 np n R y x x ... x 4 n4 n1 n2 np n R y x x ... x 5 n5 n1 n2 np the i−th subject, which is a (J × 1)-dimensional vector of all responses, where y is the binary measurement of the j-th response variable observed ij on the i-th subject. It is not difﬁcult to generalize the methodology to the case where the number of response variables differs over subjects, but that complicates the notation. As before, let x represent the multiple predictors observed on i−th subject. In Table 1, we display the structure of multivariate data in long format. The ﬁrst column, SID, is a variable which contains the subjects’ identiﬁcation number. The second column, Index, is a categor- ical indicator variable that indicates for which particular response variable the binary measurement y is obtained. In Table 1 ﬁve response variables ij are assumed, i.e., R , R ,... , R . The other columns represent measure- 1 2 5 ments of the response variable and predictor variables, respectively. A unidimensional space was used to represent the logistic regression model (3), which positions both the subjects and the two categories of the response variable. In the case of multiple responses y , the distance model can be extended to accommodate the additional responses. Suppose there is a second response variable. One possibility for generalization is to add the two points representing the categories of the second response variable to the unidimensional space. In that case, the predictor variables have a similar inﬂuence on the two response variables. 130 H.M. Worku and M. de Rooij Another generalization is that the second response variable pertains to another dimension, giving rise to a two-dimensional model. The deﬁnition of the distance becomes δ(η , γ )= (η − γ ) , im kj,m i kj m=1 where η is the coordinate for the point representing subject i on the m-th im dimension and is deﬁned as a linear combination of the predictors, η = im β + x β ;and γ is the coordinate for category k (k = {0, 1})of 0m kj,m i m response variable j on dimension m. Each response variable belongs to one and only one dimension. This assumption is driven by theories often developed by applied scientists. In the Introduction section, we discussed three different theories about comorbidity of mental disorders. Spinhoven et al. (2013), for example, found two dimensions of which the ﬁrst dimension (distress) was represented by major depressive disorder, generalized anxiety disorder, and dysthimia; and the second dimension (fear) was represented by panic disorder and social phobia. The probability for category 1 on response variable j given the pre- dictors, i.e. Pr(y =1|x )= π (x ), is now deﬁned by ij i j i exp −0.5δ(η , γ ) i 1j π (x )= . (7) j i exp −0.5δ(η , γ ) +exp −0.5δ(η , γ ) i 0j i 1j The log-odds representation of the multivariate distance model becomes, π (x ) j i 2 2 log = β (γ − γ )+ 0.5(γ − γ ) 0m 1j,m 0j,m 0j,m 1j,m 1 − π (x ) j i (8) m=1 + x β (γ − γ ) . 1j,m 0j,m i m Because each response variable belongs to a single dimension, the log odds representation can be further simpliﬁed. Suppose response variable j be- longs to dimension 1 so that γ and γ equal zero for all m> 1,i.e. 0j,m 1j,m all dimensions except the ﬁrst one. In that case, (8) simpliﬁes to a single equation instead of a sum over dimensions. This distance model for multivariate binary data (7 - 8) is called the Multivariate Logistic Distance (MLD) model. Because the MLD model is a type of bilinear model, for each dimension we have to ﬁx the origin and scale. Like in the simple logistic regression representation we ﬁx the class coordinates for one of the response variables on every dimension, e.g., γ =1 and γ =0. 1j,m 0j,m 131 Multivariate Logistic Distance Models The effect of a predictor variable on a speciﬁc response variable j is determined by the dimension the j-th response variable is positioned on. More speciﬁcally, the effect β (γ − γ ). Therefore, for dif- m 1j,m 0j,m ferent response variables on the same dimension the size of the effect is different, depending on (γ − γ ), but the direction is the same as 1j,m 0j,m long as γ ≥ γ , ∀j, m, and deﬁned by β . Furthermore, the larger 1j,m 0j,m m (γ − γ ) the bigger the effect is and vice versa. In other words, the 1j,m 0j,m larger the distance between the two points representing the categories of a single response variable, the better the predictor variables can discriminate between the categories. 2.3 Parameter Estimation The parameters in the MLD model are estimated by maximizing the likelihood function for independent data, in the multivariate situation known as quasi-likelihood; i.e., N J y (1−y ) ij ij L(θ)= π (x ) [1 − π (x )] , (9) j i j i i=1 j=1 where θ is the concatenation of all the class points and all the regression weights. Liang and Zeger (1986) showed that maximizing this quasi-likelihood provides consistent parameter estimates for the multivariate model. How- ever, the standard errors based on the corresponding Hessian matrix are bi- ased. The same authors proposed a sandwich estimator for the covariance matrix to correct for the bias (Liang and Zeger, 1986). Another method for obtaining correct standard errors is to apply a clustered bootstrap method (Sherman and Le Cessie, 1997; De Rooij and Worku, 2012; Cheng et al., 2013). In this case, the re-sampling procedure is applied on the subject (cluster) level so that the correlation between the multivariate responses is retained in each bootstrap sample. The number of independent parameters estimated in the MLD model, q, equals q =[M × (p +1)] +[(J − M) × 2]. (10) The ﬁrst term in (10), i.e., [M × (p +1)], corresponds to the number of re- gression coefﬁcients; the other term to the number of estimable class points. The identiﬁability constraints are accounted for in the second term, i.e., in each dimension the class coordinates for a single response variable are set to ﬁxed values. 132 H.M. Worku and M. de Rooij The MLD model can be ﬁtted using the NLMIXED procedure in SAS software (SAS Institute Inc., 2011). Scripts for the analyses in this paper are available upon request from the ﬁrst author. 2.4 The Relationship of the MLD Model to Generalized Estimating Equations By setting the distance between the two categories of every response variable to be equal to one, i.e., (γ − γ )=1, the MLD model 1j,m 0j,m becomes equivalent to a marginal model for multivariate binary data esti- mated using GEE method (Liang and Zeger, 1986). The restriction of the class points implies that predictor variables discriminate equally well for all response variables belonging to a speciﬁc dimension. Existing statistical packages with a GEE procedure (e.g., the genmod procedure from SAS or the geepack package from R) can be used for ﬁtting this “restricted” MLD model on multivariate binary data. Fitting the restricted MLD model using a GEE procedure involves a three-step approach: (1) construction of a design matrix for both the re- sponse and the predictor variables; and (2) applying the GEE method with the constructed design matrix; and (3) transforming the GEE parameters to MLD parameters. We now show construction of the design matrix using the example presented in Table 1. Suppose we want to ﬁt a 2-dimensional model on the ﬁve binary re- sponse variables. Further, suppose we like the ﬁrst three response variables to be represented on the ﬁrst dimension, and the fourth and the ﬁfth on the second dimension. Therefore deﬁne a response indicator matrix, denoted by (J × M). The response indicator matrix speciﬁes for Z,with dimension each response variable to which dimension it pertains, with position (j, m) equal to one if the j-th response variable belongs to the m-th dimension and zero otherwise. For the structure layed-out above, ⎡ ⎤ ⎢ 10⎥ ⎢ ⎥ Z = ⎢ 10⎥ . (11) ⎣ ⎦ The design matrix for subject i is then obtained by computing the Kronecker product between the response indicator matrix and the predictors vector (without intercept), U = Z ⊗ x , such that ⎡ ⎤ x x ... x 00 ... 0 i1 i2 ip ⎢ x x ... x 00 ... 0 ⎥ i1 i2 ip ⎢ ⎥ x x ... x 00 ... 0 U = ⎢ ⎥ . (12) i i1 i2 ip ⎣ ⎦ 00 ... 0 x x ... x i1 i2 ip 00 ... 0 x x ... x i1 i2 ip 133 Multivariate Logistic Distance Models We concatenate U and the identity matrix to get the ﬁnal design matrix, S =[I , U ], i i i ⎡ ⎤ 10 00 0 x x ... x 00 ... 0 i1 i2 ip ⎢ 01 00 0 x x ... x 00 ... 0 ⎥ i1 i2 ip ⎢ ⎥ S = 00 10 0 x x ... x 00 ... 0 . ⎢ i1 i2 ip ⎥ ⎣ ⎦ 00 01 00 0 ... 0 x x ... x i1 i2 ip 00 00 10 0 ... 0 x x ... x i1 i2 ip Then, a vertical concatenation of all S matrices will give us the ﬁnal design matrix S on which the GEE method is ﬁnally applied to obtain parameter estimates of the marginal model. This results in ﬁve response speciﬁc in- ∗ ∗ tercepts (β ,... ,β ) corresponding to the ﬁrst ﬁve columns of S and two 01 05 ∗ ∗ ∗ ∗ sets of p regression weights (β ,... ,β and β ,... ,β ). The MLD pa- 11 p1 12 p2 rameters can be derived from these as follows γ = −(β +0.5) for the 0j,m 0j dimension, m, to which disorder j belongs, zero otherwise. The regression weights β are equal to the regression weights obtained from GEE method, jm β = β . The number of parameters in the “restricted” MLD model then jm jm becomes q =[M × (p +1)] +(J − M) since additional constraints are imposed on the class points. 2.5 Model Selection In statistical analysis, we often select a parsimonious and best ﬁtting model from a set of candidate models given the data. In the MLD model, we select not only predictor variables for the ﬁnal model, but also the dimen- sionality of the model must be determined. Pan (2001) proposed the quasi-likelihood under the independence model criterion (QIC) as an extension of Akaike Information Criterion (AIC) to GEE: −1 ˆ ˆ QIC = −2L(θ)+2 trace(Ω V ), (13) where V stands for robust variance estimator obtained under the assump- tion of a general “working” covariance structure R;and Ω is for the naive variance estimator obtained under the assumption of an independence corre- lation structure. Pan (2001) also proposed a simpliﬁed version of QIC when −1 ˆ ˆ trace(Ω V ) ≈ trace(I)= q, i.e., QIC = −2L(θ)+2q. Yu and De Rooij (2013) studied the performance of QIC for de- termining the dimensionality of the Trend Vector Model (TVM). Both the 134 H.M. Worku and M. de Rooij Trend Vector model and the MLD model are marginal models in a distance framework, where the ﬁrst is used for longitudinal multinomial response variables and the latter for multivariate binary responses. Yu and De Rooij (2013) recommended QIC for determining the dimensionality of the dis- tance model. In the MLD model, we use QIC ﬁt statistics both for determining the dimensionality of the model and for variable selection. The model with the lowest QIC statistics is considered the most parsimonious and best ﬁtting model. As recommended in Yu and De Rooij (2013), we ﬁrst determine the dimensionality of the model and then proceed to the variable selection. 2.6 Biplot for the Multivariate Logistic Distance Model To enhance interpretation of the model, the results of a MLD model can be graphically represented in a biplot (Gabriel, 1971; Gower and Hand, 1996; Gower et al., 2011). The biplot represents the subjects, the response variables, and the predictor variables so that the relationship between pre- dictors and responses can be read from the graph. We ﬁrst demonstrate how the response variables are included in the biplot, and then the predictors. For a 2-dimensional MLD model the coordinates for a response vari- able are given by γ γ 0j,1 0j,2 γ = . γ γ 1j,1 1j,2 Because each response is positioned on one and only one dimension, one of the columns in γ equals zero. So, γ represents two points either on the j j ﬁrst or second dimension. Halfway between the two points, a decision line is drawn indicating equal probabilities for the two categories of a response variable. Due to these lines (horizontal for response variables on the second dimension and vertical for response variables on the ﬁrst dimension), the two dimensional space is partitioned into rectangles, each representing a most probable response proﬁle. The predictors are included in the biplot by variable axes (Gower and Hand, 1996). To derive the variable axis, ﬁrst, a pseudo-design matrix X is constructed containing ones in the ﬁrst column and zeros in the other columns except for the column representing the variable to be plotted. In this column, marker values are included within the range of the observed variable. Second, the matrix B with regression weights is deﬁned, i.e. β β 01 02 B = . β β 1 2 Finally we can compute the matrix H as H = XB, 135 Multivariate Logistic Distance Models deﬁning a straight line in our biplot. We will include variable axes for every statistically signiﬁcant predictor. Positions of the subjects are computed as the linear combination of predictor variables and are included in the biplot as points. 3. Application: The NESDA Data To illustrate the MLD model, the NESDA data (Penninx et al., 2008) introduced earlier were analysed. The sample comprised of N =2, 938 sub- jects aged 18 to 65 years (Mean = 42; S.D. =13.1). About 66.5% were female and the average number of years of education attained was 12.2 with S.D. =3.3. In this study, 37.1% of the subjects have major depressive disor- der (MDD), 10.2% have dysthmia (DYST), 15.3% have generalized anxiety disorder (GAD), 22.4% have social anxiety disorder (SP), and 28.6% have panic disorder (PD). These ﬁve disorders are the response variables. The predictors are Neuroticism (N), Extraversion (E), Openness to ex- perience (O), Agreeableness (A), and Conscientiousness (C). We also took into account three background variables, i.e., age (AGE), years of education attained (EDU), and gender (GEN: 1= female; 0= male). The linear predictor part of the MLD model is η = β + β (AGE) + β (EDU) + β (GEN) im 0m 1m i 2m i 3m i + β N + β E + β O + β A + β C , 4m i 5m i 6m i 7m i 8m i where η is a coordinate for the i-th subject position on the m-th dimen- im sion; and the β’s are regression weights. The candidate MLD models ﬁtted on the NESDA data are (1) “distress-fear” (d/f) dimensions, in which MDD, GAD, are DYST are presumed to be indicators of distress, and PD and SP for fear; (2) “depression-anxiety” (d/a) dimensions, in which MDD and DYST are indicators of depression, and GAD, PD, and SP for anxiety; (3) “unidimensional” where all the ﬁve mental disorders are indicators of a single dimension. These three theories are then translated into the following indicator matrices: ⎡ ⎤ ⎡ ⎤ ⎡ ⎤ 10 10 1 ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ 10 10 1 ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ (1) (2) (3) ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ Z = 10 , Z = 01 , Z = 1 , (14) ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ ⎣ ⎦ ⎣ ⎦ 01 01 1 01 01 1 respectively. The superscript corresponds to a theory. 136 H.M. Worku and M. de Rooij Table 2. Results of ﬁtting different MLD models to NESDA data. In the ﬁrst block, dimen- sionality of the MLD model is assessed, and followed by variable selection in the second block. Model Dimension Predictors q QIC Model Selection for Dimensionality 12(d/f) All 21 12, 396.42 22(d/a) All 21 12, 398.08 31 All 13 12, 418.87 Model Selection for Predictors 1a 2 (d/f) All 21 12396.42 1b 2 (d/f) AGE,EDU,GEN,N,E 15 12396.68 1c 2 (d/f) AGE,EDU,GEN 11 14789.41 d/f: distress/fear model. d/a: depression/anxiety model. For illustration, we ﬁtted both the MLD model with and without im- posing equal distance restrictions on the class points. The results of the MLD model with the restrictions will be presented ﬁrst, thereafter the solu- tion without the restrictions will be discussed. Table 2 shows the ﬁt statistics of the candidate MLD models. As shown in the ﬁrst block of Table 2 which compares different dimensional- ity, the 2-dimensional distress-fear (d/f) model ﬁtted the data best (QIC = 12, 396.42). In the second block of Table 2, ﬁt statistics for the compari- son of different sets of predictor variables are given. The model with all predictor variables ﬁtted the data best (QIC =12, 396.42). The regression weights of this selected model are given in Table 3. The standard errors based on both the sandwich and the clustered bootstrap method are included in Table 3. Both methods resulted in similar estimates. There is a strong positive association between neuroticism and the ˆ ˆ two dimensions: β =0.1167 with distress; and β =0.1039 with fear. 41 42 With every unit increase in neuroticism the log odds for MDD, DYST, and GAD go up by 0.1167 while the log odds for SP and PD go up by 0.1039. There is a moderate negative association between extraversion and the two ˆ ˆ dimensions: β = −0.0419 with distress; and β = −0.0320 with fear. 51 52 With every unit increase in extraversion the log odds for MDD, DYST, and GAD go down by 0.0419 while the log odds for SP and PD go down by 0.0320. From the background variables, only education has a statistically signiﬁcant effect on both dimensions: β = −0.0386 with distress; and β = −0.0575 with fear. Less educated people have a higher risk of getting a mental disorder. The variable conscientiousness had a signiﬁcant effect only on the second dimension (distress), β =0.0189, i.e. it only inﬂuences PD and SP. 137 Multivariate Logistic Distance Models Table 3. Summarized results of the ﬁnal “distress-fear” MLD model ﬁtted on NESDA data. Restriction was applied on the class points, and thus it is a restricted MLD model. The reported standard errors are based on both sandwich and clustered bootstrap methods. The number of bootstraps, B = 1000. Bootstrap Effect Parameter Estimate SE (sandwich) SE Wald Distress dimension Education β −0.0386 0.012 0.012 10.06 Gender β −0.1346 0.081 0.081 2.79 Age β 0.0012 0.003 0.003 0.15 Neuroticism β 0.1167 0.006 0.006 413.84 Extraversion β −0.0419 0.007 0.007 39.43 Openness to Experience β −0.0031 0.007 0.008 0.17 Agreeableness β −0.0074 0.008 0.007 1.03 Conscientiousness β −0.0071 0.007 0.007 1.06 Fear dimension Education β −0.0575 0.012 0.011 26.18 Gender β 0.0229 0.082 0.083 0.08 Age β −0.0008 0.003 0.003 0.08 Neuroticism β 0.1039 0.006 0.006 335.26 Extraversion β −0.0320 0.007 0.006 25.56 Openness to Experience β 0.0008 0.008 0.008 0.01 Agreeableness β −0.0003 0.008 0.008 0.00 Conscientiousness β 0.0189 0.007 0.007 6.72 statistically signiﬁcant effect, p< 0.05. Although the total number of parameters of the ﬁnal d/f model is q = 21, only sixteen of the parameters were displayed in Table 3. The others are the intercepts obtained from GEE method which are response-speciﬁc, MDD GAD DYST PD i.e., β = −2.23, β = −3.73, β = −4.28, β = −3.74,and 01 02 03 04 SP ∗ β = −4.14.Using γ = −(β +0.5) as shown in Section 2.4 and 0j,m 05 0j γ =1 + γ , the class point coordinates for each response variable 1j,m 0j,m can be obtained. Thus, γ =1.73 for MDD, γ =3.23 for GAD, 01,1 02,1 γ =3.78 for DYST, γ =3.24 for PD, and γ =3.64 for SP. 03,1 04,2 05,2 We can use the estimated class points to compare the effect of predictors on the risk of developing disorders belonging to the same dimension. For example, MDD, DYST and GAD belong to the ﬁrst dimension. Because γ =3.78 for DYST is larger than both γ =1.73 for MDD and 03,1 01,1 γ =3.23 for GAD, it means that starting from a very low subject position 02,1 on the ﬁrst dimension and then increasing this position will ﬁrst lead to higher probabilities of MDD, followed by GAD, and than for DYST. The model accounts for comorbidity because a high probability of DYST implies a highprobabilityofGAD andMDD. 138 H.M. Worku and M. de Rooij Figure 1. Biplot of the ﬁnal “distress-fear” model ﬁtted on the NESDA data. The plot is based on restrictions applied on the class points. The results of the selected MLD model are displayed in a biplot shown in Figure 1. In order to interpret the biplot, let us ﬁrst discuss how the bi- plot was constructed. The biplot is composed of two parts, i.e., the response space and the variable axes, as shown in Figures 2 and 3, respectively. The positions of the two categories of all response variables are displayed in Figure 2. For example, on the vertical dimension there are four points cor- responding to no PD, no SP, having PD, and having SP from the bottom to the top, respectively. Included in the same representation are decision lines (vertical and horizontal lines) crossing the mid-points between the two categories. The decision lines partition the two-dimensional space into rect- angles (regions), each representing a most probable response proﬁle. Each region shows the disorder proﬁle by 1’s and 0’s for the order MDD, GAD, DYST, PD, SP. An index ‘10011’, for example, corresponds to the presence of MDD, PD, and SP, but the absence of GAD and DYST. In the top-right, an index of ‘11111’ is used to represent a co-morbidity of all ﬁve mental disorders, while the region ‘00000’ in the bottom left representing the absence of disorders. In Figure 3, both the variable axes (lines) and the subjects points (grey dots) are displayed. The space includes only statistically signiﬁcant predic- 139 Multivariate Logistic Distance Models Figure 2. Representation of the binary response variables in the Euclidean space. Figure 3. Variable axes representation of the predictor variables in the Euclidean space. 140 H.M. Worku and M. de Rooij tors. On the variable axes markers are placed that represent μ ± tσ ,where x x μ is the mean of x, σ is the standard deviation, and t =0, 1, 2, 3.Variable x x labels are included at the positive side of the variable axis. Let us now interpret the biplot displayed in Figure 1. Most of the sub- jects are in the bottom left region representing absence of all the disorders. However, signiﬁcant number of subjects are in other regions representing co-morbidity of mental disorders. The regions are ‘10000’ corresponding to the presence of having only MDD; and ‘10010’ corresponding to the pres- ence of having both MDD and PD; ‘10011’ corresponding to the presence of having MDD, PD, and SP; and, ‘11011’ corresponding to the presence of all disorders, except DYST. Also a few subjects are in the upper right region having all the mental disorders. Now let us interpret the variable axes. The variable axis for Neuroti- cism (N) runs from the lower left (low values of neuroticism) to the upper right (high values of Neuroticism), indicating that persons with low values of Neuroticism are located in the ‘00000’ region, whereas persons with very high values of neuroticism are located in the ‘11111’ region. In short, the higher neuroticism the more disorders. Contrarily, the variable axes of ex- traversion points to the other direction. The length of the variable axis indicates effect size; the longer the variable axis the larger the effect of the corresponding variable on the disor- ders. The angle between the variable axis and the dimension measures the strength of their association. The smaller the angle between them, the stronger the association. For example, the angle of the extraversion vari- able axis with the ﬁrst (horizontal) dimension is relatively small compared to the angle of extraversion with the second dimension. This indicates that extraversion has a larger effect on the disorders represented on the ﬁrst di- mension (MDD, DYST, and GAD) than on the disorders presented on the second dimension (PD and SP). The angle of neuroticism with both dimen- sions is about equal, although a bit smaller with the ﬁrst dimension, indicat- ing that the relationship of neuroticism with the disorders on the ﬁrst dimen- sion (MDD, GAD, and DYST) is slightly stronger than with the other two disorders. The variable conscientiousness is highly correlated to the second dimension and not to the ﬁrst as its variable axis is orthogonal to the ﬁrst dimension. Finally, there is a strong correlation between the estimates of the sub- ject points in the two dimensions, corr(ˆ η , ηˆ )=0.99, indicating that the i1 i2 distress and fear dimensions are strongly correlated. We now present the results of MLD model that does not impose re- striction on the class points, i.e., the “unrestricted” MLD model, to address speciﬁcally the extra information from this model. The regression esti- 141 Multivariate Logistic Distance Models Table 4. Regression weights of the ﬁnal unrestricted “distress-fear” MLD model ﬁtted on NESDA data. The number of bootstraps used to obtain standard errors equals 1000. Bootstrap Effect Parameter Estimate SE Wald Distress dimension Education β −0.0203 0.006 11.45 Gender β −0.0685 0.042 2.66 Age β 0.0004 0.001 0.16 Neuroticism β 0.0605 0.004 228.77 Extraversion β −0.0226 0.004 31.92 Openness to Experience β −0.0020 0.004 0.25 Agreeableness β −0.0037 0.004 0.86 Conscientiousness β −0.0041 0.004 1.05 Fear dimension Education β −0.0202 0.005 16.32 Gender β 0.0005 0.033 0.00 Age β −0.0007 0.001 0.49 Neuroticism β 0.0424 0.003 199.75 Extraversion β −0.0141 0.003 22.09 Openness to Experience β 0.0000 0.003 0.00 Agreeableness β 0.0003 0.003 0.01 Conscientiousness β 0.0067 0.003 4.99 statistically signiﬁcant effect, p< 0.05. mates are shown in Table 4. The estimates obtained from the “unrestricted” MLD model are slightly different compared to results obtained from the “re- stricted” MLD model ﬁtted on NESDA data (shown in Table 3). However, both results lead to the same conclusion about signiﬁcance of predictors, which is also indicated by the “Wald” statistics displayed in the last column of both tables. The class points for MDD are ﬁxed for identiﬁcation on the ﬁrst dimension, i.e. the coordinates are 0 for no MDD and 1 for MDD. Similarly, the coordinates of PD on the second dimension are ﬁxed to 0 for absence and 1 for presence of the disorder. The other parameters are the class points, i.e., γ =0.96 and γ =1.73 for GAD; γ =1.10 and 02,1 12,1 03,1 γ =1.99 for DYST; and, γ = −0.25 and γ =1.28 for SP. The 13,1 05,2 15,2 distance between the two category points is 0.77 for GAD, 0.89 for DYST, and1.53for SP. This unrestricted MLD model provides additional information about how well the predictors can discriminate between the response categories. According to equation (8), the effect of the predictor variables on each re- sponse is partially determined by the distance between class points of the response variable. The larger the distance between the class points of a re- sponse variable, the better the predictor variables are able to discriminates between the categories. In the ﬁtted model, both DYST and GAD are posi- 142 H.M. Worku and M. de Rooij tioned on the ﬁrst dimension; because the distance for DYST (0.89) is larger than the distance for GAD (0.77), the effect of the predictor variables on DYST is stronger. 4. Conclusion and Discussion We proposed a multivariate logistic distance (MLD) model for ana- lyzing multivariate binary data that extends existing marginal models in a distance framework. The distance model for a single response variable was extended to analyzing multivariate binary data in the presence of predic- tors. The advantage of the MLD model over existing marginal model for multivariate data, is the possibility for dimension reduction as a form of reg- ularization which simpliﬁes the complexity of standard multivariate GLM model because less parameters are estimated. Moreover, using this dimen- sion reduction substantial theories can be represented and investigated. We have shown applications of both the “restricted” and the “unre- stricted” MLD models using an empirical data set. The former MLD model imposes a restriction on the class points and the latter model does not. The “restricted” MLD model is equivalent to a marginal model for multivari- ate binary data estimated using GEE method, which is an advantage for our model because existing software package developed for GEE can be adopted to ﬁt the MLD model. For the unrestricted case, the MLD model is a general model and can be ﬁtted by its own right. The general MLD model provides us with additional information about how well the predictors can discrimi- nate between the categories of the response variable. The MLD model has a clear interpretation where both the odds ratio expressions as well as the biplot representation can be used. The space in the biplot is partitioned into different regions that indicate the most probable response proﬁle. It is important to note that the assumption of which re- sponse variables belong to which dimension has a crucial impact on which regions might occur. In a unidimensional model there are maximal 6 re- gions, whereas in the two dimensional solution in Figure 2 there are 12 regions. Having 5 response variables, a total of 32 different proﬁles can be deﬁned. In a ﬁve dimensional model all these 32 proﬁles are present. Dimension reduction thus reduces the number of most probable response proﬁles. Moreover, the regions also account for comorbidity. In the solution of Figure 2 there is never a response proﬁle where MDD is absent and DYST and GAD are present. Similarly, if PD is present then also SP is present in the response proﬁle. Notice, however, that the model is a probability model not a deterministic model. So, a response proﬁle is most probable but the model does not say that in that region only a proﬁle must occur. The effect size of predictor variables can be read from the biplot by the length of the variable axis. The longer the variable axis the stronger the 143 Multivariate Logistic Distance Models effect. The differential effect on the two dimensions can be read from the an- gle of a variable axis with the dimension. The smaller the angle the stronger the effect. If a variable has a 90 angle with a dimension, the variable has no effect on the disorders belonging to that dimension. The MLD model is related to Canonical Correspondence Analysis (CCA), as proposed by Ter Braak (1986), which is a multivariate method used for ordination axes that maximizes the separation among the multivari- ate binary responses (Ter Braak, 1986; Ter Braak and Verdonschot, 1995). There are two main differences between CCA and our model. The ﬁrst is that these models work in different framework, i.e., the MLD model in a logistic framework where as CCA in a Gaussian framework. Due to this dif- ference, the MLD can provide a clear interpretation in terms of (log)-odds and probabilities. The second is that unlike in CCA, the MLD model can position responses (e.g., mental disorders) on certain dimensions driven by the theories that we would like to test. In areas like psychology, epidemiology, criminology, economics, po- litical sciences, etc, researchers often use Structural Equation Models (SEM) for the analysis of data similar to the NESDA data (Plewis, 1996; Von Oertzen et al., 2010; Spinhoven et al., 2013). Despite its popularity, SEM has limitations as it makes unveriﬁable assumptions about the underlying distributions of latent as well as observed variables. Moreover, SEM of- ten suffers from improper solutions, non-convergent solutions, and the pre- dicted factors are not determinate, i.e., for the same number of response variables multiple solutions can be obtained for the underlying latent vari- ables. Therefore, they cannot be uniquely identiﬁed (Acito and Anderson, 1986; Boomsma and Hoogland, 2001; Hubbard et al., 2010). In the appli- cation section, we showed that the MLD model can be used for comparing theories of interest, without making unveriﬁable assumptions about under- lying distributions. Asar and Ilk (2013) proposed marginal model with shared-parameter within the GEE method (Asar and Ilk, 2013). To compare with our MLD model, they use the ﬁve dimensional model where each response variable pertains to a unique dimension. Then, they incorporate equality restrictions for certain predictors over different dimensions, giving a so-called shared parameter. In the restricted MLD model the regression weights are shared for all response variables belonging to a speciﬁc dimension. Although our focus was on binary data, the model can be extended to polytomous data. Where for binary data there are two class points on each dimension for polytomous data there are multiple class points. Interpreta- tion follows largely the binary model, although in the ordinal case we can derive odds ratios for every contrast of two categories of a response variable. These are formed by the difference of class points coordinates, just like in 144 H.M. Worku and M. de Rooij the binary case. The polytomous situation, however, is often more compli- cated than the binary one. The binary model for every response variable is by deﬁnition unidimensional, which is not the case for polytomous data. Therefore, the polytomous case needs further study. We developed a package (an alpha version) in R,the mldm package, for ﬁtting the MLD model on multivariate binary data in the presence of predictors. The package handles both the clustered bootstrap method and the sandwich estimators for correcting standard errors of model parameters. The package provides a biplot function for the ﬁtted model. We also have SAS scripts for ﬁtting the models. The ﬁrst author can provide the package or the script upon request. References ACITO, F., and ANDERSON, R.D. (1986), “A Simulation Study of Factor Score Indetermi- nacy”, Journal of Marketing Research, 23, 111–118. AGRESTI, A. (2002), Categorical Data Analysis (2nd ed.), New York: John Wiley and Sons. AKAIKE, H. (1973), “Information Theory and an Extension of the Maximum Likelihood Principle”, in Proceedings of the Second International Symposium on Information Theory, eds. B.N. Petrov and F. Csaki, Budapest: Akademiai Kiado, pp. 267–281. ¨ ¨ ASAR, O., and ILK, O. (2013), “mmm: An R Package for Analyzing Multivariate Longitu- dinal Data with Multivariate Marginal Models”, Computer Methods and Programs in Biomedicine, 112, 649–654. BEESDO-BAUM, K. et al. (2009), “The Structure of Common Mental Disorders: A Replica- tion Study in a Community Sample of Adolescents and Young Adults”, International Journal of Methods in Psychiatric Research, 18, 204–220. BOOMSMA, A., and HOOGLAND, J.J. (2001), “The Robustness of LISREL Modeling Re- visted”, in Structural Equation Modeling: Present and Future, eds. R. Cudeck, S. de Toit, and D.Sor ¨ bom, Chicago: Scientiﬁc Software International, pp. 139–168. BORG, I. , and GROENEN, P.J.F. (2005), Modern Multidimensional Scaling: Theory and Applications (2nd ed.), New York: Springer. BULL, S.B. (1998), “Regression Models for Multiple Outcomes in Large Epidemiological Studies”, Statistics in Medicine, 17, 2179–2197. CHENG, G., YU, Z., and HUANG, J.Z. (2013), “The Cluster Bootstrap Consistency in Generalized Estimating Equations”, Journal of Multivariate Analysis, 115, 33–47. COSTA, P.T., and MCCRAE, R.R. (1992), Revised NEO Personality Inventory (NEO-PR- I) and NEO Five-Factor Inventory (NEO- FFI) Professional Manual, Odessa, FL: Psychological Assessment Resources. DE ROOIJ, M. (2009), “Ideal Point Discriminant Analysis with a Special Emphasis on Visualization”, Psychometrika, 74, 317–330. DE ROOIJ, M., and WORKU, H.M. (2012), “A Warning Concerning the Estimation of Multi- nomial Logistic Models with Correlated Responses in SAS”, Computer Methods and Programs in Biomedicine, 107(2), 341–346. ELLIOT, D.S., HUIZINGA, D., and MENARD, S. (1989), Multiple Problem Youth: Delin- quency, Substance Use, and Mental Health Problems, New York: Springer-Verlag. FITZMAURICE, G., DAVIDIAN, M., VERBEKE, G., and MOLENBERGHS, G. (2008), Longitudinal Data Analysis, London: Chapman and Hall. 145 Multivariate Logistic Distance Models GABRIEL, K.R. (1971), “The Biplot Graphical Display of Matrices with Application to Principal Component Analysis”, Biometrika, 58, 453–467. GIFI, A. (1990), Nonlinear Multivariate Analysis, Chichester: John Wiley and Sons. GOWER, J.C., and HAND, D.J. (1996), Biplots, London: Chapman and Hall. GOWER, J.C., LUBBE, S., and LE ROUX, N. (2011), Understanding Biplots, Chichester: John Wiley and Sons Ltd. HALEKOH, U., HOJSGAARD, S., and YAN, J. (2006), “The R Package geepack for Gen- eralized Estimating Equations”, Journal of Statistical Software, 15(2), 1–11. HUBBARD, A.E. et al. (2010), “To GEE or Not to GEE: Comparing Population Averaged and Mixed Models for Estimating the Associations Between Neighborhood Risk Fac- tors and Health”, Epidemiology, 21(4), 467–474. KRUEGER, R.F. (1999), “The Structure of Common Mental Disorders”, Archives of General Psychiatry, 56, 921–926. KRUSKAL, J.B., and WISH, M. (1978), Multidimensional Scaling, Sage Publications. LIANG, K.Y., and ZEGER, S.L. (1986), “Longitudinal Data Analysis Using Generalised Linear Models”, Biometrika, 73, 13–22. LIANG, K.Y., ZEGER, S.L., and QAQISH, B. (1992), “Multivariate Regression Analyses for Categorical Data, Journal of the Royal Statistical Society, Series B (Methodological), 54(1), 3–40. LIPSITZ, S.R., KIM, K., and ZHAO, L.P. (1994), “Analysis of Repeated Categorical Data Using Generalized Estimating Equations”, Statistics in Medicine, 14, 1149–1163. MCCULLAGH, P., and NELDER, J.A. (1989), Generalized Linear Models, London: Chap- man and Hall. PAN, W. (2001), “Akaike’s Information Criterion in Generalized Estimating Equations”, Biometrics, 57, 120–125. PARK, T. (1994), “Multivariate Regression Models for Discrete and Continuous Repeated Measurements”, Communications in Statistics - Theory and Methods, 23, 1547– PENNINX, B.W. et al. (2008), “The Netherlands Study of Depression and Anxiety (NESDA): Rationale, Objectives and Methods”, International Journal of Methods in Psychiatric Research, 17, 121–140. PLEWIS, I. (1996), “Statistical Methods for Understanding Cognitive Growth: A Review, A Synthesis and An Application”, British Journal of Mathematical and Statistical Psychology, 49, 25–42. R DEVELOPMENT CORE TEAM (2013), “R: A Language and Environment for Sta- tistical Computing”, Computer Software Manual Version 3.0.2, Vienna, Austria, http://www.r-project.org/. SAS INSTITUTE INC. (2011), “SAS/STAT Software”, Computer Software Manual Ver- sion 9.3, Cary, NC, http://www.sas.com. SHERMAN, M., and LE CESSIE, S. (1997), “A Comparison Between Bootstrap Meth- ods and Generalized Estimating Equations for Correlated Outcomes in Generalized Linear Models”, Communications in Statistics - Simulation and Computation, 26, 901–925. SOMMER, A., KATZ, J., and TARWOTJO, I. (1984), “Increased Risk of Respiratory Disease and Diarrhea in Children with Preexsting Mild Vitamin A Deﬁciency”, American Society for Clinical Nutrition, 40, 1090–1095. SPINHOVEN, P., DE ROOIJ, M., HEISER, W., PENNINX, B.W.J.H., and SMIT, J. (2009), “The Role of Personality in Comorbidity Among Anxiety and Depressive Disorders in Primary Care and Speciality Care: A Cross-Sectional Analysis”, General Hospital Psychiatry, 31, 470–477. 146 H.M. Worku and M. de Rooij SPINHOVEN, P., PENELO, E., DE ROOIJ, M., PENNINX, B,W., and ORMEL, J. (2013), “Reciprocal Effects of Stable and Temporary Components of Neuroticism and Affec- tive Disorders: Results of a Longitudinal Cohort Study”, Psychological Medicine, 44, 337–348. TER BRAAK, C.J.F. (1986), “Canonical Correspondence Analysis: A New Eigenvector Technique for Multivariate Direct Gradient Analysis”, Ecology, 67(5), 1167–1179. TER BRAAK, C.J.F., and VERDONSCHOT, P.F.M. (1995), “Canonical Correspondence Analysis and Related Multivariate Methods in Aquatic Ecology”, Aquatic Sciences, 57(3), 1015–1621. VAN DER HEIJDEN, P.G.M., MOOIJAART, A., and TAKANE, Y. (1994), “Correspon- dence Analysis and Contingency Models”, in Correspondence Analysis in the Social Sciences, eds. M.J. Greenacre and J. Blasius, New York: Academic Press, pp. 79– VON OERTZEN, T., HERTZOG, C., LINDENBERGER, U., and GHISLETTA, P. (2010), “The Effect of Multiple Indicators on the Power to Detect Inter-Individual Differences in Change”, British Journal of Mathematical and Statistical Psychology, 63, 627– WEI, L., and STRAM, D. (1988), “Analysing Repeated Measurements with Possibly Missing Observations by Modeling Marginal Distributions”, Statistics in Medicine, 7, 139– WEI, X. (2012), “%PROC R: A SAS Macro That Enables Native R Programming in the Base SAS Environment”, Journal of Statistical Software, 46. WORKU, H.M., and DE ROOIJ, M. (2016), “Properties of Ideal Point Classiﬁcation Models for Bivariate Binary Data”, Psychometrika (accepted for publication). ZIEGLER, A. (2011), Generalized Estimating Equations, New York: Springer. ZIEGLER, A., and ARMINGER, G. (1995), “Analyzing the Employment Status with Panel Data from GSOEP - A Comparison of the MECOSA and the GEE1 Approach for Marginal Models”, Vierteljahreshefte zur Wirtschaftsforschung, 64, 72–80. ZIEGLER, A., KASTNER, C., and BLETTNER, M. (1998), “The Generalized Estimating Equations: An Annotated Bibliography”, Biometrical Journal, 40(2), 115–139. Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org /licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.
Journal of Classification – Springer Journals
Published: Mar 26, 2018
You can share this free article with as many people as you like with the url below! We hope you enjoy this feature!
Read and print from thousands of top scholarly journals.
Already have an account? Log in
Bookmark this article. You can see your Bookmarks on your DeepDyve Library.
To save an article, log in first, or sign up for a DeepDyve account if you don’t already have one.
Copy and paste the desired citation format or use the link below to download a file formatted for EndNote
Access the full text.
Sign up today, get DeepDyve free for 14 days.
All DeepDyve websites use cookies to improve your online experience. They were placed on your computer when you launched this website. You can change your cookie settings through your browser.