Iterated filtering with a "prior"?

[walsmanjason]

I would like to use the efficiency of iterated filtering to fit a model to some field data that I have. I have laboratory data on the value of one of the fitted parameters, call it A. I know both the particle filter and iterated filtering do not support priors, which is fine, but I would like to include this information in my fitting procedure.

In terms of the code, I can easily include this information in my iterated filtering procedure by including the value of A as a data point and determining the probability of obtaining the laboratory data based on the model's estimate of the "true value" of A. There is not a true time during the data time series to include this single point of data on A so I suppose I would include it at the beginning of the time series. I understand that iterated filtering mutates parameter values over the time series so the time at which I include that A data point may matter.

My main question is whether I am committing a statistical sin that I don't see here which a reviewer may ding me for?

[kingaa]

Thanks for the interesting question, @walsmanjason. An alternative to the method you describe is to include the parameter A as a state variable, with the prior incorporated into the rinit basic component. In pomp, this is a simple matter of moving A from the parameters list to the state variables and including generation of A in the rinit code. You may of course need parameters for the prior distribution of A (which, as you no doubt know, are often called "hyperparameters"): these can simply be included in the ordinary parameter list. For the reasons you intuit in your second paragraph, this option may be more attractive than the one you propose.

As for statistical errors, one can never be sure what a reviewer will say, so my policy is to think instead about what my reply to criticism would be. I find it useful in this context to contemplate the duality between priors in the Bayesian sense and random effects in the frequentist sense. In this light, my suggestion above is to treat A as a random effect for the purposes of frequentist inference based on IF. This is mathematically equivalent to putting a prior on A.

[walsmanjason]

Thank you for the prompt reply! This makes sense (and is a good policy about reviewers!). The remaining part of my question is something that I accidentally left unspecified in my original question (mea culpa): I would like to get some information about how the value of A may be different in the field data than in the laboratory data by both incorporating my lab-based knowledge of A and fitting A to the field data. Hence my desire to basically use a prior so as to get something of a posterior distribution for A. So I see three main options for doing this:

1. My plan in the original post and be on watch for any odd behavior in the convergence of logLik or other parameters.

2. The rinit method that you describe with the value of parameter A being directly evaluated against the prior density in my dmeasure (this seems very similar to 1, I think, except that ).

3. The rinit method you describe without alterations. Then after iterated filtering gets me close to some sort of MLE, employing pmcmc with a formal prior for A (similar to what is suggested in the pomp publication). In my mind, 3 is still preferable to doing iterated filtering without any of my laboratory-based information on A; this is because I have some reason to think, from preliminary fitting, that iterated filtering would move unreasonably far away from the lab-based information on A if it does not have any constraints based on that information.

Do you have any further advice or warnings about these options? Thank you very much again!

[kingaa]

If I may summarize, it seems you have information about the parameter A from two distinct sources—the lab and the field—and you want to understand how these two perspectives on A differ. Moreover, to the extent the numerical parameter estimates do differ, it is an indication of how the situation in the field may differ from that in the lab and/or how the model may misrepresent reality. Because these issues are quite general and usually of central scientific importance, your question is particularly worthwhile.

There are two broad approaches to dealing with multiple sources of information. In the usual Bayesian approach, one wishes to assimilate all the data to obtain a posterior that summarizes the total uncertainty. The posterior then combines the data from both sources and therefore obscures discrepancies between them. To examine conflict in the data sources, one compares prior to posterior.

The second approach involves estimating the parameters separately and comparing the results directly. In effect, the model acts to translate the field data into the same terms as the lab data, facilitating comparison. In other words, one can view parameter estimates as model predictions, which become falsifiable in the presence of additional data. Let us call this second approach "Popperian".

Now both of these approaches can be accommodated by pomp. For example, in the Bayesian mode, one could implement the prior according to my suggestion, estimating the hyperparameters of that prior. If the maximum likelihood estimates of those parameters then moved far from their starting values (as set according to the lab data), this would be evidence of discrepancy between the lab and the field perspectives on the parameter A. By the way, in Bayesian parlance, the MLE of the hyperparameters obtained in this way would be called the "maximum marginal likelihood estimate".

To implement a Popperian approach, one might simply estimate the parameter A using the field data and again using the lab data (with a different model of course), and compare the two estimates. As you can see, the differences between these two approaches is not necessarily very big at the level of the algorithms.

As for the proposals you make above, let's revisit those here.

My plan in the original post and be on watch for any odd behavior in the convergence of logLik or other parameters.

What kind of "odd behavior" would one be looking for, and in comparison to what? It is far from clear that the discrepancies you are on the lookout for would manifest themselves in convergence behavior. Without some reason to think they would, this idea sounds vague and implausible.

The rinit method that you describe with the value of parameter A being directly evaluated against the prior density in my dmeasure (this seems very similar to 1, I think, except that ).

(...it seems you inadvertently dropped a phrase in your parenthesis). This proposal is unsound. Having implemented the uncertainty in the rinit component, one wouldn't want to "double count" it using the dmeasure component. That is, as we discussed before, one could either implement the uncertainty of A in rinit or in dmeasure, but not in both.

The rinit method you describe without alterations. Then after iterated filtering gets me close to some sort of MLE, employing pmcmc with a formal prior for A (similar to what is suggested in the pomp publication).

You are proposing a distinct approach here. As you say, you could do formal Bayesian inference using pmcmc, including a prior on A, which in that case could be implemented using rprior and/or dprior basic components rather than by rinit. This is distinct from the approach of using rinit to implement the random effect on A and using IF to estimate the parameters. Either would work in principle, though I suspect the latter would be more efficient.

In my mind, 3 is still preferable to doing iterated filtering without any of my laboratory-based information on A; this is because I have some reason to think, from preliminary fitting, that iterated filtering would move unreasonably far away from the lab-based information on A if it does not have any constraints based on that information.

Note that this preliminary result is evidence for the discrepancy noted above. In such cases, I always find myself wondering whether I really want to sweep the discrepancy under the rug using a fully Bayesian approach or look at it more closely to understand its source. Of course, your course of action must be guided by your purpose: is it worth trying to understand why and how the model explains the data as it does, or are you merely trying to make forecasts that are as accurate as possible?

[walsmanjason]

Thank you again for your thoughts! I should clarify and expand on some of my hasty comments above.

The overall goal of the project is to fit a small number of pomp models to field data from our wildlife host, representing different immunological hypotheses with qualitatively different effects on the distribution of the data (specifically, qualitatively different effects on the trend in the standard deviation of one of the data types). Each competing model contains the A parameter. Primarily, I want to use our lab data on A to support the model competition process by helping ensure the fitted values are in a biologically reasonable range (e.g., I happen to know the value of A in the field must be within + or - 50% of the value in the lab). Secondarily, I would like to gain some information about the difference between the value of A in the lab and the value in the field.

Regarding my suggestion of including A as a data point evaluated at the beginning of the time series in dmeasure, I did not mean to imply that the lab-field discrepancy would be revealed by the convergence behavior of log likelihood. I only meant that the artificial inclusion of A in dmeasure at a single time point in the time series may be "ok" to me as an investigator if iterated filtering still climbs up to an approximate plateau in log likelihood. Of course, your comment about implausibility may still apply if this approach makes my fitting procedure, likelihood estimate, or model competition invalid.

Circling back to your rinit suggestion, perhaps this is the way for me to go. Given hyperparameters labA, alphaA, and betaA, I could write the rinit component of A as
A=(rbeta(alphaA,betA)-0.5)*labA;

to ensure that A remains within a reasonable range. Then I would allow alphaA and betaA to change through iterated filtering.