Constructing a joint distribution from pairwise marginals

Question

Consider a set of random variables $\{X_i\}$ with joint pdf $f(x_1 ... x_n)$. Given the marginal pdfs $f_i(x_i)$, we can construct a joint distribution $$g(x_1 ... x_n) = \prod_i f_i(x_i)$$ which has the same marginals as $f$. (In particular, $g$ is the maximum entropy distribution satisfying this constraint.) $g$ is in some sense an approximation to $f$.

If we know the pairwise joint marginals $f_{ij}(x_i,x_j)$, is it possible to do something similar? Specifically:

1. What conditions do the $f_{ij}$ need to satisfy in order for there to exist a $g$ giving pairwise joint marginals $f_{ij}$?
2. Is there an explicit formula or simple algorithm that gives any such $g$?
3. Is there an explicit formula that gives the maximum-entropy such $g$?

I'm also interested in generalizations to $k$-wise joint marginals, and the conditions under which the corresponding sequence of approximations converges to $f$, but that seems like a separate question. An answer for the special case where all the $f_{ij}$ are the same would also be interesting.

See also:
Constructing a joint distribution from pairwise bivariate marginal distributions? (No answers, but a comment mentions "pair-copula model")

Answer 1

For discrete random variables the raking/iterative proportional fitting algorithm constructs a joint distribution if one exists (under some additional assumption about zero cells). It works for marginal distributions of any order, and not necessarily the same order for each margin.

IPF is (or was) used to fit loglinear models -- its guaranteed convergence rate is only linear, but it use less memory than Newton-type methods. Raking is used in survey statistics to match sample margins to known population margins.

I don't know if this generalises to continuous distributions.

Answer 2

I had the same question, and was surprised to learn the problem of combining multivariate marginal distributions to a full, joint distribution doesn't have a straightforward analytical solution, as it does with univariate marginals (assuming independence).

I don't know if iterative proportional fitting is the only or best solution, but it well worked for me. I'm posting some R code using the (unusually excellent) mipfp package with an example of implementation for posterity.

# simulate some survey data
# ethnicity (hispanic / not), sex (male / female), & race (white / other)
simdat <- array(rpois(n = 8, 10), dim= c(2, 2, 2), 
                dimnames = list(ethn = c("H", "NH"),
                                sex = c("M", "F"),
                                race = c("W", "O")))

# proportions table; compare to p.hat below
simpro <- simdat / sum(simdat) 

# all bivariate marginal distributions
sex_race <- apply(simdat, MARGIN = c(2, 3), sum)
ethn_race <- apply(simdat, MARGIN = c(1, 3), sum)
ethn_sex <- apply(simdat, MARGIN = c(1, 2), sum)
target_marginals <- list(sex_race, 
                         ethn_race, 
                         ethn_sex)

# set a seed (initial values) to all 1s 
# or use known / prior cell counts when available
seed <- (simdat + 1) / (simdat + 1)

# indicate to which indices bivariate marginals refer
# (see apply functions above. order matters!)
target_dims <- list(c(2, 3), c(1, 3), c(1, 2))

# run IPF algorithm
outp<- Ipfp(seed, target_dims, target_marginals,
               iter=50, print=TRUE, tol=1e-5)

outp$x.hat # estimated contingency table (array)
outp$p.hat # estimated multinomial probability distribution

I believe that by setting the initial values to 1, we are essentially assuming whatever independence is left over after specifying the bivariate marginals. Please correct me if I'm wrong!

Note that the algorithm converges in one iteration if you use simdat as the initial values, as expected!

Answer 3

$$ \newcommand{\avg}{\mathbb{E}} $$ The answers above are helpful for discrete distributions, but don't have much to say about continuous distributions. They also don't address the question of existence, which turns out to be pretty interesting. Here are a few things I have figured out:

1. The pairwise marginals need to result in consistent univariate distributions, e.g. we must have $$ \int_y f_{XY}(x,y) = \int_z f_{XZ}(x,z) $$
   In addition, the covariance matrix of any set of functions of the variables, which has elements given by $$ M_{ij}^{(h_1 ... h_n)} = \avg[h_i(X_i)h_j(X_j)] $$ must be positive definite. Notice that $ M_{ij} $ can be computed entirely from the pairwise marginals. I suspect that these conditions are sufficient, but haven't found a way to prove it yet. Also, I think for the second condition it might suffice to check only the functions in any complete basis set, e.g. $h_i(x) = x^{n_i}$. This case looks a little like the Hamburger Moment Problem, and it's possible that a similar proof would be useful here.
2. Under some regularity conditions, specifying all pairwise marginals is equivalent to specifying all two-variable moments $$ \mathbb{E}[x_i^n x_j^m] $$ The problem is then to make some choice for the remaining moments. This choice needs to satisfy positive-definiteness (and possibly other constraints I've missed?). This paper gives some conditions for the existence of a solution. Isserlis' theorem itself does this for the special case where the $ f_{ij} $ are Gaussian. It seems plausible that there's a generalization of Isserlis' theorem that will work for the general case, where we make some choice like like $$ \avg[X^2 Y^2 Z] = \frac{1}{3}\bigg(\avg[X^2 Y^2]\avg[Z] + \avg[X^2 Z]\avg[Y^2] + \avg[Y^2 Z]\avg[X^2]\bigg) $$
3. The maximum-entropy solution is of the form $$ g(x_1 ... x_n) = \prod_{i,j} g_{ij}(x_i,x_j) $$ and any solution of this form is maximum-entropy. This can be shown using Lagrange multipliers, and is equivalent to the log-linear model observation above. Note that Isserlis' theorem does give the maximum-entropy solution in the Gaussian case.

Answer 4

The general framework is to select from the set of all joint distributions matching a set of statistics the unique distribution with maximum entropy.

For pair-wise covariance statistics, see:
Dempster, Arthur P. Covariance selection. Biometrics (1972): 157-175.

For arbitrary "filters" (functions of many variables), see:
Zhu, S. C., Y. Wu, and D. Mumford. FRAME: Filters, Random Fields and Maximum Entropy - Towards a Unified Theory for Texture Modeling. IEEE CVPR 1993.

(Note similarities in approach. Dempster, Mumford, Wu and Zhu were all at Harvard in the early 1990s.)