A Generalized Multinomial Distribution from Dependent Categorical Random Variables

TL;DR

This paper introduces a mathematical framework for modeling dependent categorical variables using a generalized multinomial distribution, addressing correlation effects in machine learning applications.

Contribution

It provides a novel construction of dependent categorical variables and derives a generalized multinomial distribution to model their count probabilities.

Findings

Mathematical construction of dependent categorical variables

Derivation of a generalized multinomial distribution

Potential impact on machine learning algorithms

Abstract

Categorical random variables are a common staple in machine learning methods and other applications across disciplines. Many times, correlation within categorical predictors exists, and has been noted to have an effect on various algorithm effectiveness, such as feature ranking and random forests. We present a mathematical construction of a sequence of identically distributed but dependent categorical random variables, and give a generalized multinomial distribution to model the probability of counts of such variables.

Figures

Tables

Table 1. Table 1: Probability mass distribution at Level 1

${\epsilon }_1$/Branch	Path	Mass at Node	Interval
1	parent $\to$ 1	$p_1$	$(0,1/K]$
2	parent $\to$ 2	$p_2$	$(1/K,2/K]$
⋮	⋮	⋮	⋮
i	parent $\to$ i	$p_i$	$(\left(i-1\right)/K,i/K]$
⋮	⋮	⋮	⋮
K	parent $\to$ K	Mass $p_K$	$(\left(K-1\right)/K,1]$

Table 2. Table 2: Probability mass distribution at Level II, Node i

${\epsilon }_2$/Branch	Path	Mass at Node	Interval
1	$i.1\to 1$	$p_i{p}_1^{-}$	$(\frac{\left(i-1\right)K}{K^2},\frac{\left(i-1\right)K+1}{K^2}]$
2	$i.1\to 2$	$p_i{p}_2^{-}$	$(\frac{\left(i-1\right)K+1}{K^2},\frac{\left(i-1\right)K+2}{K^2}]$
⋮	⋮	⋮	⋮
i	$i.1\to i$	$p_i{p}_i^{+}$	$(\frac{\left(i-1\right)K+\left(i-1\right)}{K^2},\frac{\left(i-1\right)K+i}{K^2}]$
⋮	⋮	⋮	⋮
K	$i.1\to K$	$p_i{p}_K^{-}$	$(\frac{iK-1}{K^2},\frac{iK}{K^2}]$

Equations

$$p^{+}:=P({\epsilon }_i=1|{\epsilon }_1=1)=p+\delta q$$

$$q^{+}:=P({\epsilon }_i=1|{\epsilon }_1=0)=p-\delta p$$

$$Cor({\epsilon }_i,{\epsilon }_j)=\{\begin{array}{ll}\delta ,& i=1\\ {\delta }^2,& i\ne j,i,j\ge 2\end{array}$$

$${\epsilon }_N=i\text{ on }\left(\frac{Kj}{K^N},\frac{Kj+i}{K^N}\right]$$

$${\epsilon }_N=K\text{on }\left(\frac{K(j+1)-1}{K^N},\frac{K(j+1)}{K^N}\right]$$

$${\epsilon }_N=i$$

$${\epsilon }_N=K$$

$$p_i\prod _{j=2}^r\left[{\left(p_i^{+}\right)}^{[{\epsilon }_j=i]}\prod _{l\ne i}{\left(p_l^{-}\right)}^{[{\epsilon }_j=l]}\right],i=1,...,K$$

$$p_i=p_i{\left(p_i^{+}+\sum _{l\ne i}{p}_i^{-}\right)}^{r-1}=p_i\sum _{{\epsilon }_2,...,{\epsilon }_r}\prod _{j=2}^r\left[{\left(p_i^{+}\right)}^{[{\epsilon }_j=i]}\prod _{l\ne i}{\left(p_l^{-}\right)}^{[{\epsilon }_j=l]}\right]$$

$$P({\epsilon }_2=i)=P\left[\bigcup _{j=1}^{K^2-1}\left(\frac{Kj}{K^2},\frac{Kj+1}{K^2}\right]\right]=p_i{p}_i^{+}+p_i^{-}\sum _{j\ne i}{p}_j=p_i$$

$${\epsilon }_r=\sum _{l=1}^{K^r}{\epsilon }_r^l\text{mathds}{1}_{\left(\frac{l-1}{K^r},\frac{l}{K^r}\right]},\text{ where }{\epsilon }_r^l=\{\begin{array}{ll}l\mathrm{mod}K,& 0\not\equiv l\mathrm{mod}K\\ K,& 0\equiv l\mathrm{mod}K\end{array}$$

$$m_r^l$$

$$p_i=P({\epsilon }_{r-1}=i)=P\left(\bigcup _{\begin{array}{c}l=1\\ i\equiv l\mathrm{mod}K\end{array}}^{K^{r-1}}\{{\epsilon }_{r-1}^l=i\}\right)=\sum _{\begin{array}{c}l=1\\ i\equiv l\mathrm{mod}K\end{array}}^{K^{r-1}}P\left({\epsilon }_{r-1}^l=i\right)=\sum _{\begin{array}{c}l=1\\ i\equiv l\mathrm{mod}K\end{array}}^{K^{r-1}}{m}_{r-1}^l$$

$$m_{r-1}^l$$

$$P({\epsilon }_r=1)=\sum _{\begin{array}{c}\xi =1\\ 1\equiv \xi \mathrm{mod}K\end{array}}^{K^r}{m}_r^{\xi }=\sum _{l=1}^K{m}_{r-1}^l{p}_1^{+}+\sum _{l=K+1}^{K^{r-1}}{m}_{r-1}^l{p}_1^{-}$$

$$p_1^{+}\left(\sum _{\xi =0}^{K-1}{m}_{r-1}^{l+\xi }\right)=m_{r-1}^l{p}_1^{+}+m_{r-1}^l\left(\sum _{j=2}^K{p}_j^{-}\right)=m_{r-1}^l$$

$$P({\epsilon }_r=1)=\sum _{\begin{array}{c}\xi =1\\ 1\equiv \xi \mathrm{mod}K\end{array}}^{K^r}{m}_r^{\xi }=\sum _{l=1}^K{m}_{r-1}^l{p}_1^{+}+\sum _{l=K+1}^{K^{r-1}}{m}_{r-1}^l{p}_1^{-}=\sum _{\begin{array}{c}l=1\\ 1\equiv l\mathrm{mod}K\end{array}}^K{m}_{r-1}^l+\sum _{l=K+1}^{K^{r-1}}{m}_{r-1}^l=\sum _{l=1}^{K^{r-1}}{m}_{r-1}^l=p_1$$

$$\text{Cov}([{\epsilon }_{\iota }=i],[{\epsilon }_{\tau }=j])=\text{mathds}E[[{\epsilon }_{\iota }=i][{\epsilon }_{\tau }=j]]-\text{mathds}E[[{\epsilon }_{\iota }=i]]\text{mathds}E[[{\epsilon }_{\tau }=j]]=P({\epsilon }_{\iota }=i,{\epsilon }_{\tau }=j)-P({\epsilon }_{\iota }=i)P({\epsilon }_{\tau }=j)$$

$$\text{Cov}([{\epsilon }_1=i],[{\epsilon }_2=i])=P({\epsilon }_1=i,{\epsilon }_2=i)-p_i^2=p_i{p}_i^{+}-p_i^2=\delta p_i(1-p_i)$$

$$\text{Cov}([{\epsilon }_1=i],[{\epsilon }_2=j])=P({\epsilon }_1=i,{\epsilon }_2=j)-p_i{p}_j=p_i{p}_j^{-}-p_i{p}_j=-\delta p_i{p}_j$$

$${\epsilon }_1=1\text{ branch contains nodes }$$

$${\epsilon }_1=2\text{ branch contains nodes }$$

$$⋮$$

$${\epsilon }_1=K\text{ branch contains nodes }$$

$$P({\epsilon }_1=1,{\epsilon }_{r-1}=i)=\{\begin{array}{ll}{p}_1{p}_1^{+},& i=1\\ p_1{p}_i^{-},& i\ne 1\end{array}$$

$$p_1{p}_1^{+}=P({\epsilon }_1=1,{\epsilon }_{r-1}=1)=\sum _{\begin{array}{c}l=1\\ 1\equiv l\mathrm{mod}K\end{array}}^{K^{r-2}}{m}_{r-1}^l$$

$$p_1{p}_i^{-}=P({\epsilon }_1=1,{\epsilon }_{r-1}=i)=\sum _{\begin{array}{c}l=1\\ i\equiv l\mathrm{mod}K\end{array}}^{K^{r-2}}{m}_{r-1}^l$$

$$m_{r-1}^l=p_1^{+}{m}_{r-1}^l+m_{r-1}^l\sum _{j=2}^K{p}_j^{-}=P({\epsilon }_1=1,{\epsilon }_r=1)+m_{r-1}^l\sum _{j=2}^K{p}_j^{-}$$

$$P({\epsilon }_1=1,{\epsilon }_r=1)$$

$$=\sum _{\begin{array}{c}l=1\\ 1\equiv l\mathrm{mod}K\end{array}}^{K^{r-2}}{m}_{r-1}^l{p}_1^{+}+\sum _{j=2}^K\sum _{\begin{array}{c}l=1\\ 1\not\equiv l\mathrm{mod}K\end{array}}^{K^{r-2}}{m}_{r-1}^l{p}_j^{-}$$

$$=p_1{p}_1^{+}{p}_1^{+}+p_1{p}_1^{+}\sum _{j=2}^K{p}_j^{-}$$

$$=p_1{p}_1^{+}$$

$$P({\epsilon }_2=i,{\epsilon }_3=j)=\{\begin{array}{ll}{p}_i{\left(p_i^{+}\right)}^2+(p_i^{-}{)}^2{\sum }_{j\ne i}{p}_j,& i=j\\ p_i{p}_i^{+}{p}_j^{-}+p_j{p}_i^{-}{p}_j^{+}+{\sum }_{\begin{array}{c}l\ne i\\ l\ne j\end{array}}{p}_l{p}_i^{-}{p}_j^{-},& i\ne j\end{array}$$

$$P({\epsilon }_2=1,{\epsilon }_{\tau -1}=1)=p_1{\left(p_1^{+}\right)}^2+(p_1^{-}{)}^2\sum _{j=2}^K{p}_j^{-}$$

Full text

1. Introduction

2. 1 Background

3. 2 Construction of Dependent Categorical Variables

4. 2.1 Example construction

5. 2.2 Properties

6. 2.2.1 Identically Distributed but Dependent

7. 2.2.2 Pairwise Cross-Covariance Matrix

8. 3 Generalized Multinomial Distribution

9. 4 Properties

10. 4.1 Marginal Distributions

11. 4.2 Moment Generating Function

12. 4.3 Moments of the Generalized Multinomial Distribution

13. 5 Generating a Sequence of Correlated Categorical Random Variables

14. 5.1 Probability Distribution of a DCRV Sequence

15. 5.2 Algorithm

16. 5.2.1 Example DCRV Sequence Generation

17. 6 Conclusion

18. Acknowledgments

Introduction

Bernoulli random variables are invaluable in statistical analysis of phenomena having binary outcomes, however, many other variables cannot be modeled by only two categories. Many topics in statistics and machine learning rely on categorical random variables, such as random forests and various clustering algorithms. [6, 7]. Many datasets exhibit correlation or dependency among predictors as well as within predictors, which can impact the model used. [6, 9]. This can result in unreliable feature ranking [9], and inaccurate random forests [6].

Some attempts to remedy these effects involve Bayesian modeling [2] and various computational and simulation methods [8]. In particular, simulation of correlated categorical variables has been discussed in the literature for some time. [1, 3, 5]. Little research has been done to create mathematical framework of correlated or dependent categorical variables and the resulting distributions of functions of such variables.

Korzeniowski [4] studied dependent Bernoulli variables, formalizing the notion of identically distributed but dependent Bernoulli variables and deriving the distribution of the sum of such dependent variables, yielding a Generalized Binomial Distribution.

In this paper, we generalize the work of Korzeniowski [4] and formalize the notion of a sequence of identically distributed but dependent categorical random variables. We then derive a Generalized Multinomial Distribution for such variables and provide some properties of said distribution. We also give an algorithm to generate a sequence of correlated categorical random variables.

1 Background

Korzeniowski defined the notion of dependence in a way we will refer to here as dependence of the first kind (FK dependence). Suppose $(\varepsilon_{1},...,\varepsilon_{N})$ is a sequence of Bernoulli random variables, and $P(\varepsilon_{1}=1)=p$. Then, for $\varepsilon_{i},i\geq 2$, we weight the probability of each binary outcome toward the outcome of $\varepsilon_{1}$, adjusting the probabilities of the remaining outcomes accordingly.

Formally, let $0\leq\delta\leq 1$, and $q=1-p$. Then define the following quantities

[TABLE]

Given the outcome $i$ of $\varepsilon_{1}$, the probability of outcome $i$ occurring in the subsequent Bernoulli variables $\varepsilon_{2},\varepsilon_{3},...\varepsilon_{n}$ is $p^{+},i=1$ or $q^{+},i=0$. The probability of the opposite outcome is then decreased to $q^{-}$ and $p^{-}$, respectively.

Figure 1 illustrates the possible outcomes of a sequence of such dependent Bernoulli variables. Korzeniowski showed that, despite this conditional dependency, $P(\varepsilon_{i}=1)=p\;\forall\;i$. That is, the sequence of Bernoulli variables is identically distributed, with correlation shown to be

[TABLE]

These identically distributed but correlated Bernoulli random variables yield a Generalized Binomial distribution with a similar form to the standard binomial distribution. In our generalization, we use the same form of FK dependence, but for categorical random variables. We will construct a sequence of identically distributed but dependent categorical variables from which we will build a generalized multinomial distribution. When the number of categories $K=2$, the distribution reverts back to the generalized binomial distribution of Korzeniowski [4]. When the sequence is fully independent, the distribution reverts back to the independent categorical model and the standard multinomial distribution, and when the sequence is independent and $K=2$, we recover the standard binomial distribution. Thus, this new distribution represents a much larger generalization than prior models.

2 Construction of Dependent Categorical Variables

Suppose each categorical variable has $K$ possible categories, so the sample space $S=\{1,...,K\}$222These integers should not be taken as ordered or sequential, but rather as character titles of categories. The construction of the correlated categorical variables is based on a probability mass distribution over $K-$adic partitions of $[0,1]$. We will follow graph terminology in our construction, as this lends a visual representation of the construction. We begin with a parent node and build a $K$-nary tree, where the end nodes are labeled by the repeating sequence $(1,...,K)$. Thus, after $N$ steps, the graph has $K^{N}$ nodes, with each node labeled $1,...,K$ repeating in the natural order and assigned injectively to the intervals $\left(0,\tfrac{1}{K^{N}}\right],\left(\tfrac{1}{K^{n}},\tfrac{2}{K^{N}}\right],...,\left(\tfrac{K^{N}-1}{K^{n}},1\right]$. Define

[TABLE]

where $i=1,...,K-1$, and $j=0,1,...,K^{N-1}-1$. An alternate expression for (2) is

[TABLE]

To each of the branches of the tree, and by transitivity to each of the $K^{N}$ partitions of $[0,1]$, we assign a probability mass to each node such that the total mass is 1 at each level of the tree in a similar manner to [4].

Let $0<p_{i}<1$, $i=1,...,K$ such that $\sum_{i=1}^{K}p_{i}=1$, and let $0\leq\delta\leq 1$ be the dependency coefficient. Define

$p_{i}^{+}:=p_{i}+\delta\sum_{l\neq i}p_{l}$, $p_{i}^{-}=p_{i}-\delta p_{i}$

for $i=1,...,K$. These probabilities satisfy two important criteria:

Lemma 1**.**

• •

$\sum_{i=1}^{K}p_{i}=p_{i}^{+}+\sum_{l\neq i}p_{l}^{-}=1$**

• •

$p_{i}p_{i}^{+}+p_{i}^{-}\sum_{l\neq i}p_{l}=p_{i}$**

Proof.

The first is obvious from the definitions of $p_{i}^{+/-}$ above. The second statement follows clearly from algebraic manipulation of the definitions: $p_{i}p_{i}^{+}+p_{i}^{-}\sum_{l\neq i}p_{l}=p_{i}^{2}+p_{i}(1-p_{i})=p_{i}$ ∎

We now give the construction in steps down each level of the $K$-nary tree.

LEVEL 1:

Parent node has mass 1, with mass split $1\cdot\prod_{i=1}^{K}p_{i}^{[\varepsilon_{1}=i]}$, where $[\cdot]$ is an Iverson bracket. This level corresponds to a sequence $\varepsilon$ of dependent categorical variables of length 1.

LEVEL 2:

Level 2 has $K$ nodes, with $K$ branches stemming from each node. This corresponds to a sequence of length 2: $\varepsilon=(\varepsilon_{1},\varepsilon_{2})$. Denote $i.1$ as node $i$ from level 1. For $i=1,...,K$,

Node $\mathbf{i.1}$ has mass $p_{i}$, with mass split $p_{i}\left(p_{i}^{+}\right)^{[\varepsilon_{2}=i]}\prod\limits_{\begin{subarray}{c}j=1,j\neq i\end{subarray}}^{K}\left(p_{j}^{-}\right)^{[\varepsilon_{2}=j]}$

In this fashion, we distribute the total mass 1 across level 2. In general, at any level $r$, there are $K$ streams of probability flow at Level $r$. For $\varepsilon_{1}=i$, the probability flow is given by

[TABLE]

We use this flow to assign mass to the $K^{r}$ intervals of $[0,1]$ at level $r$ in the same way as above. We may also verify via algebraic manipulation that

[TABLE]

where the first equality is due to the first statement of Lemma 1, and the second is due to (4).

2.1 Example construction

For an illustration, refer to Figure 2. In this example, we construct a sequence of length $N=3$ of categorical variables with $K=3$ categories. At each level $r$, there are $3^{r}$ nodes corresponding to $3^{r}$ partitions of the interval $[0,1]$. Note that each time the node splits into 3 children, the sum of the split probabilities is 1. Despite the outcome of the previous random variable, the next one always has three possibilities. The sample space of categorical variable sequences of length 3 has $3^{3}=27$ possibilities. Some example sequences are $(1,3,2)$ with probability $p_{1}p_{3}^{-}p_{2}^{-}$, $(2,1,2)$ with probability $p_{2}p_{1}^{-}p_{2}^{+}$, and $(3,1,1)$ with probability $p_{3}p_{1}^{-}p_{1}^{-}$. These probabilities can be determined by tracing down the tree in Figure 2.

2.2 Properties

2.2.1 Identically Distributed but Dependent

We now show the most important property of this class of sequences– that they remain identically distributed despite losing independence.

Lemma 2**.**

$P(\varepsilon_{r}=i)=p_{i};i=1,...,K,r\in\mathbb{N}$.

Proof.

The proof proceeds via induction. $r=1$ is clear by definition, so we illustrate an additional base case. For $r=2$, from Lemma 1, and for $i=1,...,K$,

[TABLE]

Then at level $r$, in keeping with the alternative expression (3), we express $\varepsilon_{r}$ in terms of the specific nodes at level $r$:

[TABLE]

Let $m_{r}^{l}$ be the probability mass for $\varepsilon_{r}^{l}$. Then, for $l=1,...,K^{r}$ and $i=1,...,K-1$,

[TABLE]

As the inductive hypothesis, assume that for $i=1,...,K-1$,

[TABLE]

Each mass $m_{r-1}^{l}$ is split into $K$ pieces in the following way

[TABLE]

which may be written as

[TABLE]

Then

[TABLE]

When $1\equiv l\mod K$,

[TABLE]

Equation 11 holds because $m_{r-1}^{l+\xi}=\frac{p_{l+\xi}^{-}}{p_{1}^{+}}$, $\xi=0,...,K-1$, and by of Lemma 1. Thus,

[TABLE]

A similar procedure for $i=2,...,K$ follows and the proof is complete. ∎

2.2.2 Pairwise Cross-Covariance Matrix

We now give the pairwise cross-covariance matrix for dependent categorical random variables.

Theorem 1** (Cross-Covariance of Dependent Categorical Random Variables).**

Denote $\Lambda^{\iota,\tau}$ be the $K\times K$ cross-covariance matrix of $\varepsilon_{\iota}$ and $\varepsilon_{\tau}$, $\iota,\tau=1,...,n$, defined as $\Lambda^{\iota,\tau}=E[(\varepsilon_{\iota}-E[\varepsilon_{\iota}])(\varepsilon_{\tau}-E[\varepsilon_{\tau}])]$. Then the entries of the matrix are given by $\Lambda^{1,\tau}_{ij}=\begin{cases}\delta p_{i}(1-p_{i}),&i=j\\ -\delta p_{i}p_{j},&i\neq j\end{cases}$, $\tau\geq 2$, and $\Lambda^{\iota,\tau}_{ij}=\begin{cases}\delta^{2}p_{i}(1-p_{i}),&i=j\\ -\delta^{2}p_{i}p_{j},&i\neq j\end{cases}$, $\tau>\iota,$ $\iota\neq 1$.

Proof.

The $ij$th entry of $\Lambda^{\iota,\tau}$ is given by

[TABLE]

Let $\iota=1$. For $j=i$, and $\tau=2$,

[TABLE]

For $j\neq i$ and $\tau=2$,

[TABLE]

For $\tau\neq 2$, it suffices to show that $P(\varepsilon_{1}=i,\varepsilon_{\tau}=j)=\begin{cases}p_{i}p_{i}^{+},&j=i\\ p_{i}p_{j}^{-},&j\neq i\end{cases}$.

Starting from $\varepsilon_{1}$, the tree is split into $K$ large branches governed by the results of $\varepsilon_{1}$. Then at level $r$, there are $K^{r}$ nodes indexed by $l$. Each of the K large branches contains $K^{r-1}$ of these nodes. That is,

[TABLE]

Then $P(\varepsilon_{1}=1,\varepsilon_{r}^{l}=i)=m_{r}^{l};\quad i\equiv l\mod K,l=1,...,K^{r-1}$. We have already shown the base case for $r=2$. So, as an inductive hypothesis, assume

[TABLE]

Then we have that for $i=1$,

[TABLE]

and for $i\neq 1$,

[TABLE]

Moving one step down the tree (still noting that $\varepsilon_{1}=1$), we have seen that the each mass $m_{r-1}^{l}$ splits as

[TABLE]

Therefore,

[TABLE]

A similar strategy shows that $P(\varepsilon_{1}=1,\varepsilon_{r}=i)=p_{1}p_{i}^{-}$ and that $P(\varepsilon_{1}=i,\varepsilon_{r}=j)=\begin{cases}p_{i}p_{i}^{+},&j=i\\ p_{i}p_{j}^{-},&j\neq i\end{cases}$.

Next, we will compute the $ij$th entry of $\Lambda^{\iota,\tau}$, $\iota>1,\tau>\iota$. For $\iota=2,\tau=3$,

[TABLE]

This can be seen by tracing the tree construction from level 2 to level 3, with an example given in Figure 2. Next, we will show that (16) holds for $\tau>3$.

First, note that $P(\varepsilon_{2}=1,\varepsilon_{\tau}^{l}=1)=m_{\tau}^{l}$ for $1\equiv l\mod K$, and $l=(\xi-1)L^{\tau-1}+j$, where $\xi=1,...,K$, $j=1,...,K^{\tau-2}$. We have shown the base case where $\tau=3$. For the inductive hypothesis, assume that

[TABLE]

Then we have that $\left(p_{1}^{-}\right)^{2}\sum_{j=2}^{K}p_{j}^{-}=\sum_{1\equiv l\mod K}m_{\tau-1}^{l}$ for the $l$ defined above. Then, at each $m_{\tau-1}^{l}$, we have the following mass splits:

[TABLE]

Then

[TABLE]

Using the same tactic as (11), we see that using (17), adding the components, and combining the terms correctly, the proof is complete for any $\tau$. The proof for $P(\varepsilon_{2}=i,\varepsilon_{\tau}=j)$ for any $i,j$ follows similarly, and the proof for $P(\varepsilon_{\iota}=i,\varepsilon_{\tau}=j)$ follows from reducing to the above proven claims.

∎

In the next section, we exploit the desirable identical distribution of the categorical sequence in order to provide a generalized multinomial distribution for the counts in each category.

3 Generalized Multinomial Distribution

From the construction in Section 2, we derive a generalized multinomial distribution in which all categorical variables are identically distributed but no longer independent.

Theorem 2** (Generalized Multinomial Distribution).**

Let $\varepsilon_{1},...,\varepsilon_{N}$ be categorical random variables with categories ${1,...,K}$, constructed as in Section 2. Let $X_{i}=\sum_{j=1}^{N}[\varepsilon_{j}=i],i=1,...,K$, where $[\cdot]$ is the Iverson bracket. Denote $\mathbf{X}=(X_{1},...,X_{K})$, with observed values $\mathbf{x}=(x_{1},...,x_{K})$. Then

[TABLE]

Proof.

For brevity, let $\varepsilon_{(-1)}=(\varepsilon_{2},...,\varepsilon_{n})$ denote the sequence of $n$ categorical random variables with the first variable removed. Conditioning on $\varepsilon_{1}$,

[TABLE]

Now, when $\varepsilon_{1}=1$, there are ${N-1\choose x_{1}-1}$ combinations of the remaining $N-1$ categorical variables $\{\varepsilon_{i}\}_{i=2}^{N}$ to reside in category 1, ${N-1-x_{1}-1\choose x_{2}}$ ways the remaining $N-1-x_{1}-1$ categorical variables can reside in category 2, and so forth. Finally, there are ${N-1-x_{1}-1\sum_{i=2}^{K-1}x_{i}\choose x_{K}}$ ways the final unallocated $\{\varepsilon_{i}\}$ can be in category $K$. Thus,

[TABLE]

Similarly, for $i=2,...,K$

[TABLE]

Summing completes the proof. ∎

4 Properties

This section details some useful properties of the Generalized Multinomial Distribution and the dependent categorical random variables.

4.1 Marginal Distributions

Theorem 3** (Univariate Marginal Distribution).**

The univariate marginal distribution of the Generalized Multinomial Distribution is the Generalized Binomial Distribution. That is,

[TABLE]

where $q=\sum_{\begin{subarray}{c}j\neq i\end{subarray}}p_{j}$, $q^{+}=q+\delta p_{i}$, and $q^{-}=q-\delta q$

Proof.

First, we claim the following: $q^{+}=q+\delta p_{i}=p_{l}^{+}+\sum_{\begin{subarray}{c}j\neq l\\ j\neq i\end{subarray}}p_{j}^{-}$, $l=2,...,K$. This may be justified via a simple manipulation of definitions:

[TABLE]

Similarly, $q^{-}=q-\delta q$. Thus, we may collapse the number of categories to 2: Category $i$, and everything else. Now, notice that for $l\neq i$, $\left(p_{l}^{+}\right)^{x_{l}-1}\prod\limits_{\begin{subarray}{c}j\neq l\\ j\leq i\end{subarray}}\left(p_{j}^{-}\right)^{x_{j}}=\left(q^{+}\right)^{N-x_{i}-1}$ for $l=1,...,K$ and $l\neq i$. Fix $k\neq i$. Then

[TABLE]

Then

[TABLE]

∎

The above theorem shows another way the generalized multinomial distribution is an extension of the generalized binomial distribution.

4.2 Moment Generating Function

Theorem 4** (Moment Generating Function).**

The moment generating function of the generalized multinomial distribution with $K$ categories is given by

[TABLE]

where $\mathbf{X}=(X_{1},...,X_{K}),\mathbf{t}=(t_{1},...,t_{K})$.

Proof.

By definition, $M_{\mathbf{X}}(\mathbf{t})=\mathds{E}\left[e^{\mathbf{t}^{T}X}\right]=\sum_{\mathbf{X}}e^{\mathbf{t}^{T}\mathbf{X}}P(\mathbf{X}=\mathbf{x})=\sum_{\mathbf{X}}e^{\mathbf{t}^{T}\mathbf{X}}\sum_{i=1}^{K}p_{i}\frac{(N-1)!}{(x_{i}-1)!\prod_{j\neq i}x_{j}!}\left(p_{i}^{+}\right)^{x_{i}-1}\prod_{j\neq i}\left(p_{j}^{-}\right)^{x_{j}}$. Let $S_{m}=\sum_{i=1}^{m}x_{i}$, and Expanding the above,

[TABLE]

Taking the first term, denoted $T_{1}$, let $y=x_{1}-1$. Then

[TABLE]

The summation of the above is simply the moment generating function of a standard multinomial distribution with probabilities $\mathbf{p}=(p_{1}^{+},p_{2}^{-},...,p_{K}^{-})$. Thus,

[TABLE]

A similar procedure follows with the remaining terms, and summing finishes the proof.

∎

4.3 Moments of the Generalized Multinomial Distribution

Using the moment generating function in the standard way, the mean vector $\mu$ and the covariance matrix $\Sigma$ may be derived.

Expected Value

The expected value of $\mathbf{X}$ is given by $\mu=n\mathbf{p}$ where $\mathbf{p}=(p_{1},...,p_{K})$

Covariance Matrix

The entries of the covariance matrix are given by

[TABLE]

Note that if $\delta=0$, the generalized multinomial distribution reduces to the standard multinomial distribution and $\Sigma$ becomes the familiar multinomial covariance matrix. The entries of the corresponding correlation matrix are given by

[TABLE]

If $\delta=1$, the variance of $X_{i}$ tends to $\infty$ with $n$. This is intuitive, as $\delta=1$ implies perfect dependence of $\varepsilon_{2},...,\varepsilon_{n}$ on the outcome of $\varepsilon_{1}$. Thus, $X_{i}$ will either be 0 or $n$, and this spread increases to $\infty$ with $n$.

5 Generating a Sequence of Correlated Categorical Random Variables

For brevity, we will take the acronym DCRV for a Dependent Categorical Random Variable. A DCRV sequence $\varepsilon=(\varepsilon_{1},...,\varepsilon_{n})$ is in itself a random variable, and thus has a probability distribution. In order to provide an algorithm to generate such a sequence, we first derive this probability distribution.

5.1 Probability Distribution of a DCRV Sequence

The probability distribution of the DCRV sequence $\varepsilon$ of length $n$ is given formally in the following theorem. The proof follows in a straightforward fashion from the construction in Section 2 and is therefore omitted.

Theorem 5** (Distribution of a DCRV Sequence).**

Let $(\Sigma,\mathcal{F},\mathds{P})=([0,1],\mathcal{B},\mu)$. Let $\varepsilon_{i}:[0,1]\to\{1,...,K\}$, $i=1,...,n$, $n\in\mathbb{N}$ be DCRVs as defined in (2). Let $\varepsilon=(\varepsilon_{1},...,\varepsilon_{n})$ denote the DCRV sequence with observed values $e=(e_{1},...,e_{n})$. Then $\mu$ has the density

[TABLE]

and

[TABLE]

where $m^{i}$ is the mass allocated to the interval $\left(\frac{i-1}{K^{n}},\frac{i}{K^{n}}\right]$ by (4) and $i$ as the lexicographic order of $e$ in the sample space $\{1,...,K\}^{n}$ given by the relation $\frac{i}{K^{n}}=\sum_{j=1}^{K^{n}}\frac{\varepsilon_{j}-1}{K^{j}}$.

5.2 Algorithm

We take a common notion of using a uniform random variable in order to generate the desired random variable $\varepsilon$. For $\varepsilon$ with distribution $F(x)=\int_{0}^{x}f(y)dy$, $f(x)$ as in (24), it is clear that $F$ is invertible with inverse $F^{-1}$. Thus, $F^{-1}(U)$ for $U\sim\text{Uniform}[0,1]$ has the same distribution as $\varepsilon$. Therefore, sampling $u$ from $U$ is equivalent to the sample $e=F^{-1}(u)$ from $\varepsilon$.

In Section 2, we associated $\varepsilon_{n}$ to the intervals given in (2)

[TABLE]

From the construction in Section 2, each sequence has a 1-1 correspondence with the interval $\left[\frac{i-1}{K^{n}},\frac{i}{K^{n}}\right)$ for a unique $i=1,...,K^{n}$. The probability of such a sequence can be found using Theorem 5:

[TABLE]

where $l$ is the length of the above interval, and $s_{i}=\sum_{j=1}^{i}m^{j}$. Therefore, we have now partitioned the interval $[0,1)$ according to the distribution of $\varepsilon$ bijectively to the $K-$nary partition of $[0,1)$ corresponding to the particular sequence. Thus, sampling $u\in[0,1)$ from a uniform distribution and finding the interval $[s_{i-1},s_{i})$ and corresponding $i$ will yield the unique DCRV sequence.

Algorithm Strategy:

Given $u\in[0,1)$ and $n\in\mathbb{N}$, find the unique interval $[s_{i-1},s_{i})$, $i=1,...,K^{n}$ that contains $u$ by “moving down the tree” and narrowing down the “search interval” until level $n$ is reached.

We provide an explicit example prior to the pseudocode to illustrate the strategy.

5.2.1 Example DCRV Sequence Generation

Suppose $K=3$, $p_{1}=p_{2}=p_{3}=1/3$, and $\delta=1/2$, and suppose we wish to generate a DCRV sequence of length 2. Figure 3 gives the probability flow and the corresponding probability intervals $[s_{i-1},s_{i})$ that partition $[0,1)$ according to Theorem 5. Now, suppose $u=\frac{3}{4}$. We now illustrate the process of moving down the above tree to generate the sequence.

The first level of the probability tree partitions the interval $[0,1)$ into three intervals given in Figure 3. $u=\frac{3}{4}$ lies in the third interval, which corresponds to $\varepsilon_{1}=3$. Thus, the first entry of $e$ is given by $e_{1}=3$. 2. 2.

Next, the search interval is reduced to the interval from the first step $[2/3,1)$. We then generate the partitions of $[2/3,1)$ by cumulatively adding $p_{3}p_{i}^{-},i=1,2,3$ to the left endpoint $2/3$. Thus, the next partition points are

• •

$2/3+(1/3)(1/6)=13/18$,

• •

$2/3+(1/3)(1/6)+(1/3)(1/6)=7/9$, and

• •

$2/3+(1/3)(1/6)+(1/3)(1/6)+(1/3)(2/3)=1$.

Yielding the subintervals of $[2/3,1)$:

• •

$[2/3,\quad 13/18)$,

• •

$[13/18,\quad 7/9)$, and

• •

$[7/9,\quad 1)$.

We now find the interval from above that contains $u$ is the second: $[13,18,\quad 7/9)$. Thus, $\varepsilon_{2}=2$.

Since we only sought a sequence of length 2, the algorithm is finished, and we have generated the sequence $e=(3,2)$. If a longer sequence is desired, we repeat step 2 until we are at level $n$.

In general, the algorithm is given below. The $\textbf{IntervalSearch}(x,p=(p_{1},p_{2},...,p_{m}))$ procedure finds the interval $i$ built from the partitions given in the vector $p$ containing $x$ via binary search. Let $\textbf{USample}(n)$ be the procedure that samples $n$ instances of a uniformly distributed random variable. Also, let $p_{i}^{\prime}=(p_{1}^{-},...,p_{i-1}^{-},p_{i}^{+},p_{i+1}^{-},...,p_{K})$ be the “altered” probability vector for the categorical variables $2,...,n$ given $\varepsilon_{1}=i$.

6 Conclusion

Categorical variables play a large role in many statistical and practical applications across disciplines. Moreover, correlations among categorical variables are common and found in many scenarios, which can cause problems with conventional assumptions. Different approaches have been taken to mitigate these effects, because a mathematical framework to define a measure of dependency in a sequence of categorical variables was not available. This paper formalized the notion of dependent categorical variables under a first-dependence scheme and proved that such a sequence is identically distributed but now dependent. With an identically distributed but dependent sequence, a generalized multinomial distribution was derived in Section 3 and important properties of this distribution were provided. An efficient algorithm to generate a sequence of dependent categorical random variables was given.

Acknowledgments

The author extends thanks to Jason Hathcock for his suggestions and review.