2
$\begingroup$

Short version of the question: The Chinese Restaurant Process defines a distribution over partitions of $[T] := \{1, ...., T\}$. What is the expected cardinality of the $t$th block, where $t \in \{1, ..., T\}$?

Long version of the question: Suppose we have a Chinese Restaurant Process with parameter $\alpha$. As background, the $CRP(\alpha)$ specifies a conditional distribution of the form:

$$p(z_t = k|z_1, ..., z_{t-1}) = \begin{cases} \frac{\alpha}{\alpha + t - 1} & k = K_{t-1} + 1\\ \frac{N_{t-1,k}}{\alpha + t - 1} & k \leq K_{t-1} \end{cases}$$

where $N_{t-1, k}$ is the integer number of customers at the $k$th table after the $t-1$ customer has been seated, and $K_{t-1}$ is the integer number of non-empty tables.

For each partition of $T$ customers, I can write a $T$-dimensional vector $N_T$ of how many customers are at each table. For instance, if $T=3$, the sample space of $N_t$ is $\mathcal{N}_T := \{(3, 0, 0), (2, 1, 0), (1, 2, 0), (1, 1, 1) \}$. Any of these partitions has some probability mass, given by the CRP:

$$P(N_T) = \frac{\alpha^{T-1} \prod_{k=1}^K N_{T, k}!}{(\alpha+1)_{T -1 \uparrow 1}}$$

where the denominator is the so-called rising factorial, defined as: $$x_{M \uparrow a} := \prod_{m=0}^{M-1}(x + ma)$$

My question is: what is the expected value of this random $T$-dimensional vector?

$$\mathbb{E}_{N_T}[N_T] = \sum_{N_T \in \mathcal{N}_T} N_T P(N_T)$$

$\endgroup$
1
$\begingroup$

There is some ambiguity as to what is table $k$, but I am going to assume that table $k$ is the $k$th table encountered in a sequential sampling as you have mentioned. That is, $$ N_k = \sum_{t=1}^T 1\{z_t = k\} \implies \mathbb E [N_k] = \sum_{t=1}^T \mathbb P(z_t = k). $$ We also know that $$ \mathbb P(z_t = k \mid z_{t-1},\dots,z_1) = \frac{\sum_{j=1}^{t-1} 1\{z_j = k\}}{\alpha + t - 1} 1\{k \le K_{t-1}\} + \frac{\alpha}{\alpha+t-1} 1\{ k = K_{t-1}+1 \}. $$ where $K_{t-1} = \max\{z_1,z_2,\dots,z_{t-1}\}$. Taking the expectation of both sides $$ (\alpha+t-1) \mathbb P(z_t = k) = \sum_{j=1}^{t-1} \mathbb P(z_j = k \le K_{t-1}) + \alpha \mathbb P(K_{t-1} = k-1) $$ For $j \le t-1$, \begin{align*} \mathbb P(z_j = k) &= \mathbb P(z_j = k \le K_{t-1}) + \mathbb P(z_t = k, K_{t-1} < k ) \\ &= \mathbb P(z_j = k \le K_{t-1}) \end{align*} seeing that if $K_{t-1} < k$ by definition $z_j < k$ for $j \le t-1$.

Letting $\beta_{tk} := \mathbb P(z_t = k)$ and $\gamma_{tk} = \mathbb P(K_t = k)$, we have shown $$ (\alpha + t - 1) \beta_{tk} = \sum_{j=1}^{t-1} \beta_{jk} + \alpha \gamma_{t-1,k-1} $$ for $t \ge 2$. By subtracting the case of $t+1$ from that of $t$, we obtain $$ \beta_{t+1,k} - \beta_{tk} = \frac{\alpha}{\alpha + t}(\gamma_{t,k-1} - \gamma_{t-1,k-1}) $$ hence \begin{align*} \beta_{r+1,k} &= \beta_{2k} + \sum_{t=2}^r \frac{\alpha}{\alpha + t}(\gamma_{t,k-1} - \gamma_{t-1,k-1}), \quad r \ge 2, \\ \beta_{2k} &= \frac{1}{\alpha+1} 1\{k=1\} + \frac{\alpha}{\alpha+1} 1\{k=2\}. \\ \beta_{1k}&= 1\{k=1\}. \end{align*}

The problem reduces to figuring out $\gamma_{tk}$. But the distribution of $K_t$ is known as the Chinese restaurant table distribution. Its probability mass function is given as $$ \gamma_{tk} = \frac{\Gamma(\alpha)}{\Gamma(t+\alpha)} |s(t,k)| \alpha^k \cdot 1\{k \le t\}. $$ where $s(t,k)$ are the Stirling numbers of the first kind. Putting the pieces together \begin{align*} \mathbb E(N_k) &= \sum_{r=1}^T \beta_{rk} \\ &= \beta_{1k} + \beta_{2k} + \sum_{r=2}^{T-1} \beta_{r+1,k} \\ &= \beta_{1k} + (T-2) \beta_{2k} + \sum_{r=2}^{T-1} \sum_{t=1}^r \frac{\alpha}{\alpha + t}(\gamma_{t,k-1} - \gamma_{t-1,k-1}). \end{align*}


Here is a piece of R code verifying this calculation:

# library(gmp) # Needed for the Stirling1 numbers
set.seed(123)
sample_crp = function(n, a) {
  z = vector("numeric")
  z[1] = 1
  for (t in 2:n) {
    zmax = max(z)
    freq = tabulate(z, zmax)
    freq[zmax+1] = a

    z[t] = sample(zmax+1, 1, T, prob = freq)    
  }
  z
}

gam_f = function(t, k, a) {
  if (k > t) return(0)
  return( (gamma(a) / gamma(a + t)) * as.numeric(abs(gmp::Stirling1(t, k))) * a^k )
}
bet_f = function(r, k, a) {
  if (r == 1) {
    return(if (k==1) 1 else 0)
  }
  if (r == 2) {
    if (k==1) return(1/(a+1))
    if (k==2) return(a/(a+1))
    return(0)
  }
  temp = bet_f(2, k, a)
  for (t in 2:(r-1)) {
    temp = temp + a*(gam_f(t, k-1, a) - gam_f(t-1, k-1, a))/(a + t)
  }
  temp
}
N_f = function(n, k, a) sum(sapply(1:n, function(r) bet_f(r, k, a)))

nrep = 2000
n = 20
a = 0.6
X = matrix(nrow = nrep, ncol = n)
for (rep in 1:nrep) {
  X[rep, ] = sample_crp(n, a)  
}

tab = apply(X, 1, FUN = function(z) tabulate(z, n))
rbind(estim=rowMeans(tab)[1:7], true=sapply(1:7, function(k) N_f(n, k, a)))

The output is

         [,1]     [,2]     [,3]      [,4]      [,5]      [,6]        [,7]
estim 12.8615 4.770500 1.694000 0.5070000 0.1355000 0.0290000 0.002500000
true  12.8750 4.772263 1.670688 0.5159738 0.1329717 0.0277614 0.004646026

Some details about the above argument:

  • What happened when you took the expectation of both sides? In general, the conditional expectation $\mathbb E(X\mid Y)$ is a random variable (a function of $Y$). As a random variable, $\mathbb E(X\mid Y)$ itself has an expectation and we have the general "smoothing" or "tower" property: $$ \mathbb E[\mathbb E(X\mid Y)] = \mathbb E(X). $$ Now take $X = 1_{A}$ to be the indicator of a set, i.e., $X(\omega) = 1_A(\omega) = 1\{\omega \in A\}$. Then, applying the above and noting that the expectation of $1_A$ is the probability of $A$, we have $$ \mathbb E [ \mathbb P(A \mid Y) ] = \mathbb P(A). $$ Note that $\mathbb P(A \mid Y)$ is a random variable and not a deterministic quantity. This type of argument is what is going on in that step.
$\endgroup$
8
  • $\begingroup$ Wonderful! Thank you! Yes, you correctly deduced what I mean by table $k$. $\endgroup$ – Rylan Schaeffer Jan 31 at 21:15
  • $\begingroup$ When you write, "Taking the expectation of both sides," your left-hand side is $\mathbb{P}(z_t)$. Could you clarify what happened there? It seems to me you switched from the conditional distribution $\mathbb{P}(z_t | z_{<t})$ to the marginal distribution $\mathbb{P}(z_t)$, but how is that accomplished by taking an expectation? $\endgroup$ – Rylan Schaeffer Jan 31 at 21:32
  • $\begingroup$ I'm also confused by something right after "We also know that". Your left hand side, the conditional distribution $\mathbb{P}(z_t∣z_{<t})$ is entirely deterministic, right? That is, given $z_{<t}$, we know exactly how much mass to place on each outcome for $z_t$. But you then say to take the expectation of both sides. How does one take the expectation of a deterministic function? Your next step suggests to me that you're treating the indicators as random variables - is this correct? If so, the left hand side is deterministic but the right hand side is random - how can this be? $\endgroup$ – Rylan Schaeffer Jan 31 at 21:47
  • 1
    $\begingroup$ @RylanSchaeffer, I have added some comments to end of the answer. Please have a look and let me know if that makes it clear. $\mathbb P (z_t = k \mid z_{<t})$ is not deterministic (in general) but itself a random variable. It is a function of $z_{< t}$ as you can also see from the RHS of the equation. ($z_{<t}$ is random.) $\endgroup$ – passerby51 Jan 31 at 23:47
  • $\begingroup$ One quick clarification question: in your expression for $\beta_{r+1, k}$, why is $\beta_{1, k}$ not included in the sum? $\endgroup$ – Rylan Schaeffer Feb 2 at 16:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.