Franks 2019 - Logistic Selection with Exponential Family Mixtures

Summary

If the log-odds of treatment are linear in a sufficient statistic $s(Y(t))$ of the potential outcome and the observed outcome is an exponential family (or a mixture of them), then the missing-outcome distribution is available in closed form: the natural parameter shifts by ${\gamma }_t$ and mixture weights are re-weighted. Condition 1 holds automatically and the normalising constant is tractable, so the observed-data model is fit once and the sensitivity analysis is a post-processing step.

Overview

Logistic selection (Eq. 7)

$$f(T=1|Y(t),X)={\text{logit}}^{-1}\{{\alpha }_t(X)+{\gamma }_t^{\prime }{s}_t(Y(t))\}.$$

$\gamma =({\gamma }_0,{\gamma }_1)$ are the sensitivity parameters (how treatment depends marginally on each potential outcome); ${\alpha }_t(X)$ is identified by the data once ${\gamma }_t$ is set. Parameters are assumed independent of covariates (extension: ${\gamma }_t(X)$).

Properties: for any ${\gamma }_t$ giving a proper missing distribution, a unique ${\alpha }_t(X)$ satisfies Condition 1 (Robins et al. 2000); the missing distribution is free of ${\alpha }_t(X)$ (Rotnitzky et al. 2001; Scharfstein et al. 2003):

$$f^{mis}(Y(t)|T=1-t,X)=f^{obs}(Y(t)|T=t,X)\frac{\exp ({\gamma }_t^{\prime }{s}_t(Y(t)))}{C({\gamma }_t,X)},C({\gamma }_t,X)=\int e^{{\gamma }_t^{\prime }{s}_t(y)}{f}^{obs}(y|T=t,X)dy(8)$$

and ${\gamma }_0,{\gamma }_1$ are variation independent. The difficulty is computing $C$; Scharfstein et al. used the empirical CDF or MCMC on the joint space; the logistic-mEF class avoids it.

Sign convention

As printed, Eqs. 8, 10, 12 and 22 use $+{\gamma }_t$ in the tilt for both arms. Derivation from (7): the odds of the other arm is $f(T=1-t|y)/f(T=t|y)$, which equals $e^{+({\alpha }_0+{\gamma }_{0}y)}$ for $t=0$ but $e^{-({\alpha }_1+{\gamma }_{1}y)}$ for $t=1$. So with selection defined as $P(T=1|Y(t))$ for both arms, the missing control mean shifts by $+{\gamma }_0{\sigma }_0^2$ and the missing treated mean by $-{\gamma }_1{\sigma }_1^2$ in the normal case. The interpretation in Sec. 5.1 (Examples 5-7: ${\gamma }_1>0$ means treated observed mean exceeds the unobserved) and the results in Fig. 2 are consistent with this strict version, so the printed ”+” for $t=1$ is best read as a notational shortcut where ${\gamma }_1$ is defined through selection into the arm. Check your implementation with a simulation (the code below does this).

Main Content

Proposition 1 (single exponential family)

Let $f_t^{obs}(Y(t)|T=t,X)=h(Y(t))g({\eta }_t(X))e^{s(Y(t){)}^{\prime }{\eta }_t(X)}$ with natural parameter ${\eta }_t(X)$. Under logistic selection (7), if ${\eta }_t^{\ast }(X)={\eta }_t(X)+{\gamma }_t$ lies in the natural parameter space, the missing outcomes are in the same family,

$$f_{{\gamma }_t,t}^{mis}(Y(t)|T=1-t,X)=h(Y(t))g({\eta }_t^{\ast }(X))e^{s(Y(t){)}^{\prime }{\eta }_t^{\ast }(X)},C({\gamma }_t,X)=\frac{g({\eta }_t(X))}{g({\gamma }_t+{\eta }_t(X))}.$$

Proposition 2 (mixtures of exponential families)

Let $f_t^{obs}={\sum }_k{\pi }_k{h}_k(Y(t))g_k({\eta }_{tk})e^{s(Y(t){)}^{\prime }{\eta }_{tk}}$ with common sufficient statistic $s$. Under logistic selection with ${\eta }_{tk}^{\ast }={\eta }_{tk}+{\gamma }_t$ inside each component’s natural parameter space, the missing distribution is a $K$-component mixture of the same families,

$$f^{mis}=\sum _k{\pi }_k^{\ast }{h}_k(Y)g_k({\eta }_{tk}^{\ast })e^{s(Y{)}^{\prime }{\eta }_{tk}^{\ast }},{\pi }_k^{\ast }=\frac{{\pi }_k{g}_k({\eta }_k)/g_k({\eta }_k^{\ast })}{{\sum }_{k^{\prime }}{\pi }_{k^{\prime }}{g}_{k^{\prime }}({\eta }_{k^{\prime }})/g_{k^{\prime }}({\eta }_{k^{\prime }}^{\ast })}(12\text{-}13),$$

and $C({\gamma }_t,X)={\sum }_k{\pi }_k{g}_k({\eta }_{tk})/g_k({\gamma }_t+{\eta }_{tk})$. Includes Dirichlet-process mixtures and BART observed models.

Worked cases

Example 2: binary outcome

$f_t^{obs}=\text{Bern}({\text{logit}}^{-1}{\mu }_t(X))$ and selection $\propto \text{Bern}({\text{logit}}^{-1}({\alpha }_t(X)+{\gamma }_{t}Y(t)))$. Then the unobserved arm is $\text{Bern}({\text{logit}}^{-1}({\mu }_t(X)+{\gamma }_t))$: an additive log-odds shift.

Example 3: normal outcome, quadratic selection

$f_t^{obs}=N({\mu }_t(X),{\sigma }_t^2)$, sufficient statistics $(Y,Y^2)$, selection ${\text{logit}}^{-1}({\alpha }_t(X)+{\gamma }_{t}Y(t)+{\psi }_{t}Y(t{)}^2)$. Then

$$f^{mis}=N\left(\frac{{\mu }_t(X)+{\gamma }_t{\sigma }_t^2}{1-2{\psi }_t{\sigma }_t^2},\frac{{\sigma }_t^2}{1-2{\psi }_t{\sigma }_t^2}\right)(11).$$

With ${\psi }_t=0$ the standard deviations of observed and missing distributions coincide (mean shift ${\gamma }_t{\sigma }_t^2$). Propriety requires ${\psi }_t<1/(2{\sigma }_t^2)$ (footnote 3).

Example 4: normal mixture

$f_t^{obs}={\sum }_k{\pi }_{k}N({\mu }_{tk}(X),{\sigma }_{tk}^2)$, ${\psi }_t=0$: $f^{mis}={\sum }_k{\pi }_k^{\ast }N({\mu }_{tk}+{\gamma }_t{\sigma }_{tk}^2,{\sigma }_{tk}^2)$ with

$${\pi }_k^{\ast }\propto {\pi }_k\exp \{{\gamma }_t{\mu }_{tk}(X)+\frac{1}{2}{\gamma }_t^2{\sigma }_{tk}^2\}={\pi }_k{E}_k\left[e^{{\gamma }_{t}Y}\right].$$

(This is the expansion of $g_k({\eta }_k)/g_k({\eta }_k^{\ast })$ for a Gaussian component, i.e. the moment generating function of component $k$ at ${\gamma }_t$. The paper prints a different-looking exponent, $\frac{1}{2}\{{\mu }_{tk}^2/{\sigma }_{tk}^2-({\mu }_{tk}/{\sigma }_{tk}-{\gamma }_t{)}^2\}$, which does not match this derivation dimensionally, so the form above is this note’s own derivation.) Both the mixture weights and the component means move with ${\gamma }_t$: for ${\gamma }_t>0$ components with larger mean or variance are up-weighted.

The same machinery handles semi-continuous data (zero-inflated income, see Franks 2019 - NHANES and Job Training Applications) and any Dirichlet-process mixture of exponential families.

Examples

import numpy as np
# Verify the Gaussian tilt: sample from the observed arm, importance-weight with exp(+-gamma*y), compare to closed form
rng = np.random.default_rng(1)
mu, sig, gamma = 3.0, 2.0, 0.15
y = rng.normal(mu, sig, 1_000_000)
for sign, name in [(+1, "control arm  (missing = treated)"), (-1, "treated arm  (missing = control)")]:
    w = np.exp(sign*gamma*y); w /= w.sum()
    print(name, "weighted mean", (w*y).sum(), " closed form", mu + sign*gamma*sig**2)

PyMC-style sketch (untested): fit the observed arm once, then for each posterior draw compute the shifted complete-data mean:

# mu_post: posterior draws of mu_t(x_i) per unit, sigma_post: draws of sigma_t, p_t: share of units observed in arm t
# ate_draw(gamma0, gamma1) = [p1*mu1_obs + (1-p1)*(mu1_miss)] - [p0*mu0_obs + (1-p0)*mu0_miss]
# mu1_miss = mu1 - gamma1*sigma1**2   (strict sign)    mu0_miss = mu0 + gamma0*sigma0**2

Connections

• Tilts correspond to partial-pooling free “bias parameters” in Rosenman 2020 - Combining Observational and Experimental Datasets - Overview, where the bias is estimated from an experiment instead of postulated.
• BART observed model: Chipman 2010 - BART - Overview; DP mixtures: Teh 2010 - Dirichlet Process - Overview.

See Also

Franks 2019 - Flexible Sensitivity Analysis Without Observable Implications - Overview, Franks 2019 - Calibrating Sensitivity Parameters