I'm working on fitting an exponential model
$\mathrm{Flux} = A+Bt+F\left(\exp(t_0-t/T_r) + \exp(t-t_0/T_f)\right)^{-1}+…$
to astronomical data (a light curve).
$A$, $B$, $F$, $t_0$, $T_r$, and $T_f$ are all free parameters in the model.
I am trying to determine if a two-peak, three-peak, or four-peak model is appropriate for my dataset (with support from statistical tests). For each n-peak model (n=2,3,4) there is a corresponding exponential function in my model with it's own values of $t_0$ and F. I'm using nls in R to perform the fit.
The problem I'm having is that it is ambiguous whether a 3-peak or 4-peak model is preferential. I also looked at a 5 peak model just for curiosity's sake. So far I have looked at a likelihood ratio test, the AIC and BIC and none seem to put preference on one model versus the other.
Likelihood ratio test
Using the 2-peak model as my null hypothesis
3 peak: p=1.325564e-10
4 peak: p=2.220446e-16
5 peak: p=1.110223e-16
Here I see it as the 4-peak model being preferred over the 2-peak model, as the $H_0$ is rejected with a greater significance. But when you get p values this small, do differences this small really matter?
3-peak model as $H_0$
4 peak: 5.164630e-08
5 peak: 1.877226e-08
My interpretation: 4-peak model preferred.
4-peak model as $H_0$
5 peak: .01661512
My interpretation: 5-peak model preferred over 4-peak model.
AIC
2 peak:-58011.87
3 peak:-58053.36
4 peak:-58082.92
5 peak:-58087.11
My interpretation: 4-peak and 5-peak models equally preferred. But I would select the 4-peak model over the 5-peak model because it has fewer free parameters. But I question this because the difference between the AIC values (even between the 2-peak and 5-peak models) is so small.
BIC
2 peak: -57988.89
3 peak: -58025.27
4 peak: -58049.72
5 peak: -58048.80
My interpretation: 4-peak model and 5-peak model equally preferred. Again, I would select the 4-peak model over the 5-peak model because it has fewer free parameters (also the BIC value is slightly smaller)
I understand that when performing a LRT on nonlinear models it is preferable to use a non-parametric fit for your null hypothesis, however I have been unable to get that to work in R. I think my dataset may be too small (I have 95 data points).
So my question really boils down to: have I interpreted these statistical tests properly.
Also, are there any methods better than the ones I've already utilized to select one model over the other?
Just for some context, the reason the model selection is so vital is that the values for $T_f$ under the 3-peak versus 4-peak model are very different and point to two differing physical explanations.
Best Answer
If you are trying to do model selection, I would definetly go for AIC and/or BIC.
The main reason for not looking at the likelihood-ratio tests that you gave in your example is that the p-values aren't comparable with each other (well...in general, p-values aren't comparable!): each time you calculate a p-value, you are not considering that there exist other possible models (e.g., when you perform the 2-peak model v/s the 3, 4, etc. separately) and therefore the null distribution that serves for the calculation of the p-value does not make sense for the problem at hand...the p-values are, therefore, meaningless for your problem!
Now, in general, to compare the AIC values in order to select a model, a difference of more than $2$ between the AIC values is usually denoted as "distinguishable models" by the AIC (see this wikipedia article on the AIC). Due to the low number of samples, however, I would recommend calculating the AICc instead of the "plain AIC", which is the "low-sample corrected" version of the AIC and is given by: $$AICc=AIC+\frac{2k(k+1)}{n-k-1}$$ where $AIC$ is the AIC value that you have, $k$ is the number of parameters and $n$ is the number of datapoints. You can check out more about the AICc in the wikipedia article that I cited before. For further references, please check out the excellent textbook of Burnham & Anderson (2002).