andrewcparnell
commented
7 years ago Hi kabrantes,
Apologies for only answering this now. Somehow I managed to miss it.
If you have very few observations then running a model with a residual error (as is the default) is dangerous as it can't estimate the between observation standard deviations. In fact when I run your code above I get a warning:
Warning message:
In simmr_mcmc(simmr_in) :
At least 1 group has less than 4 observations - either put each observation in an individual group or use informative prior informationYou should follow the advice in the warning! Hope that helps,
Andrew
Hi Andrew, I've got yet another question about simmr. For example, if I run this model (script below) with 2 mixture samples, I get a contibution of seagrass of 3-81% (95% CI). This lower bound of the 95% CI sounds a bit low to me given the position of the mixtures on the d13C-d15N plane... Anyway, If I then go and run the model using simmr solo based on 1 only sample, that has d13C and d15N that correspond to the average of those 2 initial samples (d13C = -11.4, d15N=8.15), I get a much more constrained result on seagrass contribution: 64-84%. Looking at the plot, that one makes more sense. I realise the source polygon in this case is a bit funny, but can't explain the large difference in results, can you help please? THANK YOU!
mix = matrix(c(-11.8, -11.0, 8.1, 8.2), ncol=2,nrow=2) colnames(mix) = c('d13C','d15N') s_names = c("Seagrass", "Mangroves","MPB", "Saltmarsh") s_means = matrix(c(-10.1,-27.3,-19.2,-25.7,1.4,2.3,3.4,6.0), ncol=2, nrow=4) s_sds = matrix(c(0.3, 1.5, 1.5, 0.7, 2.4, 0.2, 1.5, 2.1), ncol=2, nrow=4) c_means = matrix(c(1.7, 1.7, 1.7, 1.7, 5.6, 5.6, 5.6, 5.6), ncol=2, nrow=4) c_sds = matrix(c(0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5), ncol=2, nrow=4)
simmr_in = simmr_load(mixtures=mix, source_names=s_names, source_means=s_means, source_sds=s_sds, correction_means=c_means, correction_sds=c_sds) simmr_out = simmr_mcmc(simmr_in)
plot(simmr_in,xlab=expression(paste(delta^13, "C (\u2030)",sep="")), ylab=expression(paste(delta^15, "N (\u2030)",sep="")), title='Aldrichetta (P Davis) Summer')
summary(simmr_out,type='statistics') summary(simmr_out,type='quantiles')