# Regression analysis and parameter estimates with populations

I've seen a little bit here about the difference between statistical inference for random samples, and what happens when we actually have population data. Most arguments seem to suggest you "never actually have the population" and the population data you think you have represents some unobservable super population which is the data-generating process.

But say, we have a state, covered by counties, and we have data on some environmental aspect of each county, for example, the total area of forest. We want to see if some binary county-level outcome is related to forest - let's say presence or absence of a disease (thinking Lyme disease here). As far as I can tell, this qualifies as "population" level data. Would running a regular logistic regression, and the normal tests on parameters, yield correct estimates?

One thing that comes into this is whether the "samples" are independent. In this case, forest in county A is probably related to forest in the county B next door, meaning a reduction in the degrees of freedom and imprecision in standard error estimates (as I understand it). But even if we had totally independent data, would standard, sample-based statistical inference be appropriate here?

I like what chl has written. Despite that, I am moved to discuss whether this situation necessarily requires a complicated model. But first, let's begin with responses to some comments.

(1) You don't lose any degrees of freedom due to spatial correlation of forest cover. This is an explanatory variable, not the one you're trying to model. You might lose "degrees of freedom" when the residuals of the dependent variable exhibit spatial autocorrelation. Such correlation is not necessarily a given even when a map of the dependent variable itself suggests strong spatial correlation. The reason is that the correlation in the map likely derives from the correlation in the forest cover (and other spatially distributed covariates). Remember, in such models you don't ask whether the "samples are independent"--it's usually obvious they are not--but whether any deviations in them from their modeled values are independent.

(2) Therefore, a conditional autoregressive model might not be necessary. I think this modeling choice would be attractive only if you want to test a theory of contagion.

Now to answer the original question: yes, first run an ordinary logistic (or Poisson) model, because as a general principle it's good to try simple (but reasonable) models first. Provided its residuals do not exhibit strong spatial correlation, you are then ok with the results. If there is evidence of correlation and evidence that it will appreciably affect your answers (coefficients, predictions, or whatever), consider using a generalized linear geostatistical model (GLM). These are described in Diggle & Riberio, Model-based Geostatistics (a relatively inexpensive and accessible text), which itself documents several R packages for making the estimates and some associated EDA tools (geoR and geoRglm). The GLM approach lets you simultaneously fit your model and assess the degree of spatial autocorrelation. The principal limitations I have found in these packages are (1) they don't handle anisotropy well--you can detect it but it's hard to incorporate it in a model--and (2) they don't have a provision for nested variograms, which somewhat limits your ability to model the spatial correlation. Both of these are no problem for smallish datasets, because you need (typically) hundreds to thousands of observations, or more, to model correlation at this level of detail.

Finally, a word about the "population" question. I assume you are interested in more than a mere description of the data: you seek information about a possible association between disease and other observable factors. Even when you have a comprehensive description of the data for a spatial region, it still does not act like a census, because the outcomes could have turned out otherwise. Next year, with identical forest cover, there will be a slightly different pattern of disease. In other regions of the country or the world and at other times, exactly the same combinations of explanatory variable values are likely to produce varying rates of disease. Thus, you're modeling a process, not a population.

• Thanks for your (always) enlightened remarks! You provide a more thorough response than mine wrt. OP (I must admit I'm a little bit biased toward epidemiology). [I'll +1 ASAP] – chl Sep 26 '10 at 20:40
• Thanks a lot for this explanation. It really helps. The distinction between process and population is what I was missing (among other stuff!).. that the process is stochastic, and will vary from place to place and time to time – kip Sep 26 '10 at 21:47

What you describe seems to refer to a particular case of multilevel modeling, where data are organized into a hierarchical structure; in your case, forests (1st level unit) nested in counties nested in states, but see Under what conditions should one use multilevel/hierarchical analysis?.

Now, the "particular" case comes from the fact that you want to account for the spatial proximity which might be a vector for the propagation of Lyme disease (correct me if I am wrong), as is done in epidemiology where one is interested in studying the geography of infectious disease. In the usual case, we can use so-called spatial models like the multiple membership model or the conditional autoregressive model, among others. I enclose a couple of references about these approaches at the end, but I think you will find more references by looking at related studies in ecology or epidemiology.

Now, I think that you may pay a particular attention at the following paper of Langford et al. which features multilevel modeling with spatially correlated data:

Langford, IH, Leyland, AHL, Rasbash, and Goldstein, H (1999). Multilevel modelling of the geographical distributions of diseases. Journal of Royal Statistical Society C, 48, 253-268.

Harvey Goldstein is the author of an excellent book on multilevel modeling, Multilevel Statistical Models (the 2nd edition is available for free). Finally, the book of Andrew Gelman, Data Analysis Using Regression and Multilevel/Hierarchical Models, may provide additional clues about hierarchical/multilevel modeling.

About software, I know there is the R spdep package for modeling spatially correlated outcomes, but there are some examples of analysis of spatial hierachical data with WinBUGS on the BUGS Project.

References

1. Browne, W.J., Goldstein, H. and Rasbash, J. (2001) Multiple membership multiple classification (MMMC) models. Statistical Modelling, 1, 103-124.
2. Lichstein, JW, Simons, TR, Shriner, SA, and Franzreb, KE (2002). Spatial autocorrelation and autoregressive models in ecology. Ecological Monographs, 72(3), 445-463.
3. Feldkircher, M (2007). A Spatial CAR Model applied to a Cross-Country Growth Regression.
4. Lawson, AB, Browne, WJ, and Vidal Rodeiro, CL (2003). Disease mapping with WinBUGS and MLwiN. John Wiley & Sons.
• Thanks for both of your responses. I will dig a little further into some of the references you suggested chi. I guess what I was getting at in my question, was aside from the issue of spatial effects, if I actually did have population-wide independent data, would classical hypothesis tests suffice - since they are based on sampling theory? Is multi-level modelling the only approach here? Forests aren't really nested in counties, rather I'd say they make up a measured proportion of each county, and those measurements are correlated in space... – kip Sep 26 '10 at 20:16
• @kip Would you mind explaining "they make up a measured proportion of each county"? I initially thought there was a % of Lyme disease recorded at each site (in your case all the forests in each country). – chl Sep 26 '10 at 20:29
• @kip About your second point, well if you consider you have a set of measurements collected on a fixed period of time on all possible statistical units (this is what you called your population), you can apply inferential procedures provided you are willing to assume there is an evolving generating model underlying the observed data (but see the related question I put in comment for a thorough discussion, e.g. by @ars, on this point). – chl Sep 26 '10 at 20:29
• @chi Sorry for the confusion, I meant a proportion for forest cover, and a binary indicator for presence or absence at each country. I think I get what you are saying here, would the assumption of not having an evolving generating model essentially be an assumption of 'stationarity'? – kip Sep 26 '10 at 21:41
• @kip Finally, I like @whuber point of view: "you're modeling a process, not a population", which seems pretty close to what I had in mind. – chl Sep 27 '10 at 9:27

I work in a related area using what is called the Population Approach, in pharmacometrics. Basically you have a sample of individuals from the general population, and a moderate-size sample of observations of each individual, plus demographic covariates for these individuals.

To model them, we want to make levels of model. We want to make an overall model that applies to the individuals that we hope will be representative of the entire population. The parameters of this overall model are called "fixed effects".

Below that level, we have a model of each individual, because we expect individuals to differ from each other. So we have additional parameters that are different for each individual, and these are called "random effects".

A key part to be estimated in the overall model is a measure of the variability of the random effects. I.e. we know the individuals will differ, but we want to model that variability.

A further level of modeling is the variability of individual observations.

Applying this to your field is not something I fully understand, but I would recommend WinBugs as a modeling tool, and a book "Markov Chain Monte Carlo in Practice" by Gilks, Richardson, and Spiegelhalter. There are example problems in there that look like they might be applicable to your problem.