Hey–here’s an R package for imputing Census data using iterative proportional fitting on available margins.

Gustavo pointed us to this project (see also here) by David Dorer to impute data to lower census levels (block/tract) from PUMAs (public use microdata areas, as defined by the U.S. Census).

From a quick glance this looks like the basic method that I would recommend right now. I have dreams of a better, Bayesian, approach, but iterative proportional fitting (IPF) is what I usually recommend to people. So if they’ve programmed this up and set it up with real census data, that looks like it could be useful.

Ben Goodrich adds:

I think it sounds sort of reasonable. In general with the PUMS, the person-weights only depend on race, sex, and age, so I don’t know how great it is for other variables anyway. But if it rakes in such a way that the smaller subareas compose into PUMAs, that is a start. I think if you did table fusion on a bunch of variables including the geographic subarea, then you could loop through each row in the PUMS can use as a “prior” the proportion of the person’s PUMA in each subarea and use as the “likelihood” the probability of having the person’s demographics given that they live in in each subarea with positive prior probability, you could get a “posterior” probability that the live in each subarea and impute from that categorical distribution. I doubt it would matter that much, except when the subareas are Congressional districts, though.

Going forward, I’d like to move to soft constraints, as discussed here and section 2 of this post. I’m thinking of a model in which the specified margins are reals, not integers, with normal errors on the margins with sd’s that are user-specified (I think we could come up with reasonable defaults, too, given that the method would often be used in default settings).

2 thoughts on “Hey–here’s an R package for imputing Census data using iterative proportional fitting on available margins.”

  1. The undefined table fusion in the quote is something I have been working on to condition on a set of related American Community Survey (ACS) tables when the data generating process is governed by a Dirichlet-multinomial distribution, which is closed under aggregation and conditioning so the log-probability of observing the values in that set of ACS tables can be derived (i.e. in Stan) from a coherent joint distribution for the union of variables in the ACS tables.

  2. I’ve been talking to Andrew about this (Yajuan Si’s continuing to work on this, too). It’s pretty easy to do what you need for poststratification by treating the contingency table entries as a simplex (what you actually need for post stratification) and then treating the margins as multinomial observations given the marginal probabilities defined from the simplex. You could generalize from multinomial to normal, but I like that the multinomial cooks in counts. If you code this up in Stan, it’s then very easy to add other soft constraints or distributional assumptions, such as row covariance or column covariance or both with Kronecker structure or more interacted. You can also fit when the margins are not from the same surveys and don’t have the same size. Here’s an example of what I’m talking about that just lets you impose a multivariate logistic prior on the entries.

    functions {
      vector row_margin(matrix x) {
        return x * rep_vector(1, cols(x));
      }
      vector column_margin(matrix x) {
        return (rep_row_vector(1, rows(x)) * x)';
      }
    }
    data {
      int M, N;
      vector[M * N] m;
      cov_matrix[M * N] S;
      array[M] int y_row;
      array[N] int y_col;
    }
    transformed data {
      matrix[M * N, M * N] L_S = cholesky_decompose(S);
    }
    parameters {
      vector[M * N] theta_unc;
    }
    transformed parameters {
      matrix[M, N] theta = to_matrix(softmax(theta_unc), M, N);
    }
    model {
      theta_unc ~ multi_normal_cholesky(m, L_S);
    
      y_row ~ multinomial(row_margin(theta));
      y_col ~ multinomial(column_margin(theta));
    }
    

    It’s not the most efficient way to compute because it reallocates 1 vectors each time for the marginals. One could also decompose the covariance matrix down into a Kronecker product of column and row covariances with matrix normal structure. The resulting products are much easier to solve since you can solve a Kronecker matrix by taking the Kronecker product of the solutions of its components (but you never ever want to actually build the huge Kronecker matrix as I’ve done in the example code). In Stan, it’s even simpler—you just collect the rows into an array and give them the column covariance and collect the columns into an array and apply the row covariance.

    I’m generally receptive to making discrete things continuous and generally find it a good idea to listen to Andrew on modeling, but in this case I’m not sure what he’s thinking in terms of a continuous generalization. I think we could just take a multivariate normal approximation to the multinomial the same way I can take a univariate normal approximation to the binomial, but I’ve never tried it, and don’t know if it’s what Andrew intended.

Leave a Reply

Your email address will not be published. Required fields are marked *