I have a couple of ordinary differential equation models that I'm trying to fit to time-dependent biological data ($y_n$). One model is more complex than then other as it has more free parameters.

I would like to test, which model is better overall for fitting the data and estimate which model is more 'likely' to be a better reflection of the dynamics (as much as models can reflect the dynamics in biological systems)

I could use a simple AIC (corrected as there's only ~25 points) comparing the two fits and calculate the Aikaike weights, but how this performs is dependent on the uncertainty in the data, as any effect could well be small and swamped by technical noise in the data.

To combat this I have generated bootstrap time-series samples using the Gaussian process-based bootstrap (as in http://www.ncbi.nlm.nih.gov/pubmed/19289448) giving ($\hat{y}^i_n$), and fitted the two models to each bootstrapped dataset. This then gives me a residual sum of squares (RSS) and an AIC(c) value for each bootstrap sample.

Can I use the distributions of the AIC(c)s as a measure of which model is more likely? The distributions are not independent, so an individual bootstrap-sample-based AIC(c) comparison could give me the information I want (scoring each model based on the number of times its AIC(c) value is the lower of the two), but I wonder if there's anything more formal.


It seems I'm answering my own questions now. Anyhow, according to Burnham & Anderson (Model selection and multimodel inference: a practical information-theoretic approach), comparing AIC (or equivalent information criterion) per bootstrapped sample and tallying up a model frequency (no of selections / no of bootstraps) per model is probably the best option.

On the data I was working on (which I haven't quoted), this does give some very obvious model choices (which back up the model choices from the original datasets), where the model frequencies for one particular model are well above half the no of bootstraps. The "chosen" model in each case differs depending on the data set used but, for our purposes, this is ideal as it allows us to classify data sets in terms of which model is the best fit overall (and make some guess as to the underlying biology).

  • 1
    $\begingroup$ For real world fits, I prefer to have at least 5 samples per parameter. I love having a few hundred, but when things start getting sparse, I consider that a line into danger. You have very few samples. I hope your equations have few parameters. $\endgroup$ – EngrStudent - Reinstate Monica Mar 26 '16 at 23:10
  • 2
    $\begingroup$ I've got 200 GPR bootstrap series per dataset (with 300 more generated, ready to fit if this proves not powerful enough), and the models have only 2,3 or 4 parameters, so the ratio is at least 50 series per parameter. $\endgroup$ – user36196 Mar 28 '16 at 23:13

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.