Dimension-Grouped Mixed Membership Models for Multivariate Categorical Data

In this subsection, we briefly review the existing MMM literature to give our proposal appropriate context. Let $K$ be the number of extreme latent profiles. Denote the $K$-dimensional probability simplex by $\Delta^{K-1}=\{(\pi_{1},\ldots,\pi_{K}):\pi_{k}\geq 0\text{ for all }k,\;\sum_{k=1}^{K}\pi_{k}=1\}$. Each subject $i$ has an individual proportion vector $\boldsymbol{\pi}_{i}=(\pi_{i,1},\ldots,\pi_{i,K})\in\Delta^{K-1}$, which indicates the degrees to which subject $i$ is a member of the $K$ extreme profiles. The general mixed membership models summarized in Airoldi et al., (2014) have the following distribution,

where $\boldsymbol{\pi}_{i}=(\pi_{i,1},\ldots,\pi_{i,K})$ follows the distribution $D_{\boldsymbol{\alpha}}$ and is integrated out; the $\boldsymbol{\alpha}$ refers to some generic population parameters depending on the specific model. The hierarchical Bayesian representation for the model in (1) can be written as follows.

where “i.i.d.” is short for “independent and identically distributed”. The number $p$ in (1) is the number of “characteristics”, and $R$ is the number of “replications” per characteristic. As shown in (1), for each characteristic $j$, there are a corresponding set of $K$ conditional distributions indexed by parameter vectors $\{\boldsymbol{\lambda}_{j,k}:\,k=1,\ldots,K\}$. Many different mixed membership models are special cases of the general setup (1). For example, the popular Latent Dirichlet Allocation (LDA) (Blei et al., , 2003; Blei, , 2012; Anandkumar et al., , 2014) for topic modeling takes a document $i$ as a subject, and assumes there is only $p=1$ distinct characteristic (one single set of $K$ topics which are distributions over the word vocabulary) with $R>1$ replications (a document $i$ contains $R$ words which are conditionally i.i.d. given $\boldsymbol{\pi}_{i}$); LDA further specifies $D_{\boldsymbol{\alpha}}(\boldsymbol{\pi}_{i})$ to be the Dirichlet distribution with parameters $\boldsymbol{\alpha}=(\alpha_{1},\ldots,\alpha_{K})$.

Focusing on MMMs for multivariate categorical data, there are generally many characteristics with $p\gg 1$ and one replication of each characteristic with $R=1$ in (1). Each variable $y_{i,j}\in\{1,\ldots,d_{j}\}$ takes one of $d_{j}$ unordered categories. For each subject $i$, the observables $\boldsymbol{y}_{i}=(y_{i,1},\ldots,y_{i,p})^{\top}$ are a vector of $p$ categorical variables. MMMs for such data are traditionally called Grade of Membership models (GoMs) (Woodbury et al., , 1978). GoMs have been extensively used in applications, such as disability survey data (Erosheva et al., , 2007), scholarly publication data (Erosheva et al., , 2004), and data disclosure risk and privacy (Manrique-Vallier and Reiter, , 2012). GoMs are also useful for psychological measurements where data are Likert scale responses to psychology survey items, and educational assessments where data are students’ correct/wrong answers to test questions (e.g. Shang et al., , 2021).

In GoMs, the conditional distribution $f(y_{i,j}\mid\boldsymbol{\lambda}_{j,k})$ in (1) can be written as $\mathbb{P}(y_{i,j}\mid\boldsymbol{\lambda}_{j,k})=\prod_{c_{j}=1}^{d_{j}}\lambda_{j,c_{j},k}^{\mathbb{I}(y_{i,j}=c_{j})}$. Hence, the probability mass function of $\boldsymbol{y}_{i}$ in a GoM is

The hierarchical Bayesian representation for the model in (2) can be written as follows.

See a graphical model representation of the GoM with sample size $n$ in Figure 1(b), where individual latent indicator variables $(z_{i,1},\ldots,z_{i,p})\in[K]^{p}$ are introduced to better describe the data generative process.

We emphasize that the case with $p>1$ and $R=1$ is fundamentally different from the topic models with $p=1$ and $R>1$, with the former typically involving many more parameters. This is because the “bag-of-words” assumption in the topic model with $R>1$ disregards word order in a document and assumes all words in a document are exchangeable. In contrast, our mixed-membership model for multivariate categorical data does not assume a subject’s responses to the $p$ items in a survey/questionnaire are exchangeable. In other words, given a subject’s mixed membership vector $\boldsymbol{\pi}_{i}$, his/her responses to the $p$ items are independent but not identically distributed (because they follow categorical distributions governed by $p$ different sets of parameters $\{\boldsymbol{\lambda}_{j,k}\in\mathbb{R}^{d}:~{}k\in[K\}$ for $j=1,\ldots,p$); whereas in a topic model, given a document’s latent topic proportion vector $\boldsymbol{\pi}_{i}$, the $p$ words in the document are independent and identically distributed, following the categorical distribution with the same set of parameters $\{\boldsymbol{\lambda}_{k}\in\mathbb{R}^{V}:~{}k\in[K]\}$ (here $V$ denotes the vocabulary size). In this sense, the GoM model has greater modeling flexibility than topic models and are more suitable for modeling item response data, where it is inappropriate to assume that the items in the survey/questionnaire are exchangeable or share the same set of parameters. This fact is made clear also in Figure 1(b), where for each $j\in[p]$ there is a population quantity, the parameter node $\boldsymbol{\Lambda}_{j,:,:}$ (also denoted by $\boldsymbol{\Lambda}_{j}$ for simplicity), that governs its distribution. Thus identifiability is a much greater challenge for GoM models. To our best knowledge, the identifiability issue of the grade-of-membership (GoM) models for item response data considered in Woodbury et al., (1978) and Erosheva et al., (2007) has not been rigorously investigated so far. Motivated by the difficulty of identifying GoM in its original setting due to the large parameter complexity, we next propose a new modeling grouping component to enhance identifiability. Our resulting model still does not make any exchangeability assumption of the items, but rather leverages the variable grouping structure to reduce model complexity.

Generalizing Grade of Membership models for multivariate categorical data, we propose a new structure that groups the $p$ observed variables in the following sense: any subject’s latent membership is the same for variables within a group but can differ across groups. To represent the key structure of how the $p$ variables are partitioned into $G$ groups, we introduce a notation of the grouping matrix $\mathbf{L}=(\ell_{j,g})$. The $\mathbf{L}$ is a $p\times G$ matrix with binary entries, with rows indexed by the $p$ variables and columns by the $G$ groups. Each row $j$ of $\mathbf{L}$ has exactly one entry of “1” indicating group assignment. In particular,

Our key specification is the following generative process in the form of a hierarchical Bayesian representation,

where “ind.” is short for “independent”, meaning that conditional on $z_{i,g}=k$, subject $i$’s observed responses to items in group $g$ are independently generated. Hence, given the population parameters $(\mathbf{L},\boldsymbol{\Lambda},\boldsymbol{\alpha})$, the probability distribution of $\boldsymbol{y}_{i}$ can be written as

For a sample with $n$ subjects, assume the observed responses $\boldsymbol{y}_{1},\ldots,\boldsymbol{y}_{n}$ are independent and identically distributed according to the above model.

We visualize the proposed model as a probabilistic graphical model to highlight connections to and differences from existing latent structure models for multivariate categorical data. In Figure 1, we show the graphical model representations of two popular latent structure models for multivariate categorical data in (a) and (b), and for the newly proposed Gro-M³ in (c) and (d). The $\boldsymbol{\Lambda}_{j}$ for $j\in[p]$ denotes a $d_{j}\times K$ matrix with entries $\lambda_{j,c_{j},k}$. Each column of $\boldsymbol{\Lambda}_{j}$ characterizes a conditional probability distribution of variable $y_{j}$ given a particular extreme latent profile. We emphasize that the variable grouping is done at the level of the latent allocation variables $z$, and that the $\mathbf{\Lambda}_{j}$ parameters are still free without constraints just as they are in traditional LCMs or GoMs. From the visualizations in Figure 1 we can also easily distinguish our proposed model from another popular family of methods, the co-clustering methods (Dhillon et al., , 2003; Govaert and Nadif, , 2013). Co-clustering usually refers to simultaneously clustering the subjects and clustering the variables, where subjects within a cluster exhibit similar behaviors and variables within a cluster also share similar characteristics. Our Gro-M³, however, does not restrict the $p$ variables to have similar characteristics within groups, but rather allows them to have entirely free parameters $\boldsymbol{\Lambda}_{1},\ldots,\boldsymbol{\Lambda}_{p}$. The “dimension-grouping” happens at the latent level by constraining the latent allocations behind the $p$ variables to be grouped into $G$ statuses. Such groupings give rise to a novel probabilistic hybrid tensor decomposition visualized in Figure 1(c)–(d); see the next Section 2.3 for details.

Other than facilitating model identifiability (see Section 3), our dimension-grouping modeling assumption is also motivated by real-world applications. In general, our new model Gro-M³ with the variable grouping component can apply to any multivariate categorical data to simultaneously model individual mixed membership and capture variable similarity. For example, Gro-M³ can be applied to survey/questionnaire response data in social sciences, where it is not only of interest to model subjects’ partial membership to several extreme latent profiles, but also of interest to identify blocks of items which share similar measurement goals within each block. We next provide numerical evidence to demonstrate the merit of the variable grouping modeling component. For a dataset simulated from Gro-M³ (in the setting as the later Table 2) and also the real-world IPIP personality test dataset (analyzed in the later Section 6), we calculate the sample Cramer’s V between item pairs. Cramer’s V is a classical measure of association between two categorical variables, which gives a value between 0 and 1, with larger values indicating stronger association. Figure 2 presents the plots of the sample Cramer’s V matrix for the simulated data and the real IPIP data. This figure shows that the pairwise item dependence for the Gro-M³-simulated data looks quite similar to the real-world personality test data. Indeed, after fitting the Gro-M³ to this IPIP personality test dataset, the estimated model-based Cramer’s V shown in Figure 2(c) nicely and more clearly recovers the item block structure. We conjecture that many real world datasets in other applied domains exhibit similar grouped dependence.

The Gro-M³ class has interesting connections to popular tensor decompositions. For a subject $i$, the observed vector $\boldsymbol{y}_{i}$ resides in a contingency table with $\prod_{j=1}^{p}d_{j}$ cells. Since the MMMs for multivariate categorical data (both traditional GoM and the newly proposed Gro-M³) induce a probability of $\boldsymbol{y}_{i}$ being in each of these cells, such probabilities $\{p\left(y_{i,1}=c_{1},\ldots,y_{i,p}=c_{p}\mid-\right);\,c_{j}\in[d_{j}]~{}\text{for each}~{}j\in[p]\}$ can be arranged as a $p$-way $d_{1}\times d_{2}\times\cdots\times d_{p}$ array. This array is a tensor with $p$ modes and we denote it by $\mathbf{P}$; Kolda and Bader, (2009) provided a review of tensors. The tensor $\mathbf{P}$ has all the entries nonnegative and they sum up to one, so we call it a probability tensor. We next describe in detail the tensor decomposition perspective of our model; such a perspective will turn out to be useful in the study of identifiability.

The probability mass function of $\boldsymbol{y}_{i}$ under the traditional GoM model can be written as follows by exchanging the order of product and summation,

Then $\mathbf{\Phi}^{\text{{GoM}}}:=\left(\phi^{\text{{GoM}}}_{k_{1},\ldots,k_{p}};\;k_{j}\in[K]\right)$ forms a tensor with $p$ modes, and each mode has dimension $K$. Further, this tensor $\mathbf{\Phi}$ is a probability tensor, because $\phi_{k_{1},\ldots,k_{p}}\geq 0$ and it is not hard to see that its entries sum up to one. Viewed from a tensor decomposition perspective, this is the popular Tucker decomposition (Tucker, , 1966); more specifically this is the nonnegative and probabilistic version of the Tucker decomposition. The $\mathbf{\Phi}^{\text{{GoM}}}$ represents the Tucker tensor core, and the product of $\{\lambda_{j,c_{j},k}\}$ form the Tucker tensor arms.

It is useful to compare our modeling assumption to that of the the Latent Class Model (LCM; Goodman, , 1974), which follows the graphical model shown in Figure 1(a). The LCM is essentially a finite mixture model assuming each subject $i$ belongs to a single cluster. The distribution of $\boldsymbol{y}_{i}$ under an LCM is

Based on the above definition, each subject $i$ has a single variable $z_{i}\in[K]$ indicating which latent class it belongs to, rather than a mixed membership proportion vector $\boldsymbol{\pi}_{i}$. Denoting $\boldsymbol{\nu}^{\text{LC}}=(\nu_{k};\,k\in[K])$, then (6) corresponds to the popular CP decomposition of tensors (Hitchcock, , 1927), where the CP rank is at most $K$.

Finally, consider our proposed Gro-M³,

where $f(y_{i,j}|\lambda_{j,c_{j},k})$ generally denotes the conditional distribution of variable $y_{i,j}$ given parameter $\lambda_{j,c_{j},k}$. In our Gro-M³, $\lambda_{j,c_{j},k}$ specifically refer to the categorical distribution parameters for the random variable $y_{i,j}$; that is, $\lambda_{j,c_{j},k}=\mathbb{P}(y_{i,j}=c_{j}\mid z_{i,j}=k)$ denotes the probability of responding $c_{j}$ to item $j$ given that the subject’s realization of the latent profile for item $j$ is the $k$th extreme latent profile. In this case, $\mathbf{\Phi}^{\text{Gro-M${}^{3}$}}:=\left(\phi^{\text{Gro-M${}^{3}$}}_{k_{1},\ldots,k_{G}};\;k_{g}\in[K]\right)$ forms a tensor with $G$ modes, and each mode has dimension $K$. There still is $\sum_{k_{1}=1}^{K}\cdots\sum_{k_{G}=1}^{K}\allowbreak\phi_{k_{1},\ldots,k_{G}}^{\text{Gro-M${}^{3}$}}\allowbreak=1$. This reduces the size of the core tensor in the classical Tucker decomposition because $G<p$. The Gro-M³ incorporates aspects of both the CP and Tucker decompositions, providing a probabilistic hybrid decomposition of probability tensors. The CP is obtained when all variables are in the same group, while the Tucker is obtained when each variable is in its own group; see Figure 1 for a clear illustration of this fact.

Gro-M³ is conceptually related to the collapsed Tucker decomposition (c-Tucker) of Johndrow et al., (2017), though they did not model mixed memberships, used a very different model for the core tensor $\boldsymbol{\Phi}$, and did not consider identifiability. Nonetheless and interestingly, our identifiability results can be applied to establish identifiability of c-Tucker decomposition (see Remark 1 in Section 4). Another work related to our dimension-grouping assumption is Anandkumar et al., (2015), which considered the case of overcomplete topic modeling with the number of topics exceeding the vocabulary size. For such scenarios, the authors proposed a “persistent topic model” which assumes the latent topic assignment persists locally through multiple words, and established identifiability. Our dimension-grouped mixed membership assumption is similar in spirit to this topic persistence assumption. However, the setting we consider here for general multivariate categorical data has the multi-characteristic single-replication nature ($p>1$ and $R=1$); as mentioned before, this is fundamentally different from topic models with a single characteristic and multiple replications ($p=1$ and $R>1$).

For LDA and closely related topic models, there is a rich literature investigating identifiability under different assumptions (Anandkumar et al., , 2012; Arora et al., , 2012; Nguyen, , 2015; Wang, , 2019). Typically, when there is only one characteristic ($p=1$), $R\geq 2$ is necessary for identifiability; see Example 2 in Wang, (2019). However, there has been limited consideration of identifiability of mixed membership models with multiple characteristics and one replication, i.e., $p>1$ and $R=1$. GoM models belong to this category, as does the proposed Gro-M³s, with GoM being a special case of Gro-M³s.

We consider the general setup in (1), where $\boldsymbol{\Phi}$ denotes the $G$-mode tensor core induced by any distribution $D(\boldsymbol{\pi}_{i})$ on the probability simplex $\Delta^{K-1}$. The following definition formally defines the concept of strict identifiability for the proposed model.

Definition 1 gives the strongest possible notion of identifiability of the considered population quantities $(\mathbf{L},\boldsymbol{\Lambda},\boldsymbol{\Phi})$ in the model. In particular, the strict identifiability notion in Definition 1 requires identification of both the continuous parameters $\boldsymbol{\Lambda}$ and $\boldsymbol{\Phi}$, and the discrete latent grouping structure of variables in $\mathbf{L}$. The following theorem proposes sufficient conditions for the strict identifiability of Gro-M³s.

Theorem 1 requires that each of the $G$ variable groups contains at least three variables, and that for each of these $3G$ variables, the corresponding conditional probability table $\boldsymbol{\Lambda}_{j}$ has linearly independent columns. Theorem 1 guarantees not only the continuous parameters are identifiable, but also the discrete variable grouping structure summarized by $\mathbf{L}$ is identifiable. This is important practically as typically the appropriate variable grouping structure is unknown, and hence needs to be inferred from the data.

The conditions in Theorem 1 essentially requires at least $3G$ conditional probability tables, each being a matrix of size $d_{j}\times K$, to have full column rank. This implicitly requires $d_{j}\geq K$. Tan and Mukherjee, (2017) proposed a moment-based estimation approach for traditional mixed membership models and briefly discussed the identifiability issue, also assuming $d_{j}\geq K$ with some full-rank requirements. However, the cases where the number of categories $d_{j}$’s are small but the number of extreme latent profiles $K$ is much larger can arise in applications; for example, the disability survey data analyzed in Erosheva et al., (2007) and Manrique-Vallier, (2014) have binary responses with $d_{1}=\cdots=d_{p}=2$ while the considered $K$ ranges from 2 to 10. Our next theoretical result establishes weaker conditions for identifiability that accommodates $d_{j}<K$, by taking advantage of the dimension-grouping property of our proposed model class.

Before stating the theorem, we first introduce two useful notions of matrix products. Denote by $\bigotimes$ the Kronecker product of matrices and by $\bigodot$ the Khatri-Rao product. Consider two matrices $\mathbf{A}=(a_{i,j})\in\mathbb{R}^{m\times r}$, $\mathbf{B}=(b_{i,j})\in\mathbb{R}^{s\times t}$; and another two matrices $\mathbf{C}=(c_{i,j})=(\boldsymbol{c}_{\boldsymbol{:},1}\mid\cdots\mid\boldsymbol{c}_{\boldsymbol{:},k})\in\mathbb{R}^{n\times k}$, $\mathbf{D}=(d_{i,j})=(\boldsymbol{d}_{\boldsymbol{:},1}\mid\cdots\mid\boldsymbol{d}_{\boldsymbol{:},k})\in\mathbb{R}^{\ell\times k}$, then there are $\mathbf{A}\bigotimes\mathbf{B}\in\mathbb{R}^{ms\times rt}$ and $\mathbf{C}\bigodot\mathbf{D}\in\mathbb{R}^{n\ell\times k}$ with

The above definitions show the Khatri-Rao product is the column-wise Kronecker product. The Khatri-Rao product of matrices plays an important role in the technical definition of the proposed dimension-grouped MMM. The following Theorem 2 exploits the grouping structure in $\mathbf{L}$ to relax the identifiability conditions in the previous Theorem 1.

Compared to Theorem 1, Theorem 2 relaxes the identifiability conditions by lifting the full-rank requirement on the individual matrices $\boldsymbol{\Lambda}_{j}$’s. Rather, as long as the Khatri-Rao product of several different $\boldsymbol{\Lambda}_{j}$’s have full column rank as specified in (10), identifiability can be guaranteed. Recall that the Khatri-Rao product of two matrices $\boldsymbol{\Lambda}_{j_{1}}$ of size $d_{j_{1}}\times K$ and $\boldsymbol{\Lambda}_{j_{2}}$ of size $d_{j_{2}}\times K$ has size $(d_{j_{1}}d_{j_{2}})\times K$. So intuitively, requiring $\boldsymbol{\Lambda}_{j_{1}}\bigodot\boldsymbol{\Lambda}_{j_{2}}$ to have full column rank $K$ is weaker than requiring each of $\boldsymbol{\Lambda}_{j_{1}}$ and $\boldsymbol{\Lambda}_{j_{2}}$ to have full column rank $K$. The following Example 2 formalizes this intuition.

Example 2 shows that the Khatri-Rao product of two matrices seems to have full rank under fairly mild conditions. This indicates that the conditions in Theorem 2 are much weaker than those in Theorem 1 by imposing the full-rankness requirement only on a certain Khatri-Rao product of the $\boldsymbol{\Lambda}_{j}$-matrices, instead of on individual $\boldsymbol{\Lambda}_{j}$s. To be more concrete, the next Example 3 illustrates Theorem 2, as a counterpart of Example 1.

Example 2 shows that the Khatri-Rao product of conditional probability tables easily has full column rank in a toy case, and Example 3 leverages this observation to establish identifiability for almost all parameters in the parameter space using Theorem 2. We next generalize this observation to derive more practical identifiability conditions, under the generic identifiability notion introduced by Allman et al., (2009). Generic identifiability generally means that the unidentifiable parameters belong to a set of Lebesgue measure zero with respect to the parameter space. Its definition adapted to the current Gro-M³s is given as follows.

Compared to the strict identifiability notion in Definition 1, the generic identifiability notion in Definition 2 is less stringent in allowing the existence of a measure zero set of parameters where identifiability does not hold; see the previous Example 2 for an instance of a measure-zero set. Such an identifiability notion usually suffices for real data analyses (Allman et al., , 2009). In the following Theorem 3, we propose simple conditions to ensure generic identifiability of Gro-M³s.