SpBayes for an offset model

I am running spBayes to fit an 'offset' model y ~ 1.

I have a dataframe like this

    ID   lon   lat      y
1   A  90.0  5.9  0.957096100
2   A  90.5  6.0  0.991374969
3   A  91.1  6.0  0.991374969
4   A  92.7  6.1  0.913501740
5   A  94.0  6.1  0.896575928
6   A  97.8  5.2  0.631320953
7   A  98.9  4.4 -0.282432556
8   A 101.2  2.8  1.829053879
9   A 102.3  2.0  0.993621826
10  A 105.8  0.5  0.038677216

where the variable ID is a factor with two levels A and B. I would like to find a offset for the two IDs. However, when I run

fit.by.ALL <- spLM(formula=y ~ ID, data= df, coords=coords, 
               priors=priors, tuning=tuning, starting=starting,
               cov.model = "exponential", n.samples = n.samples, verbose = TRUE, 
               n.report = 50)

which gives the result

Iterations = 1:251
Thinning interval = 1 
Number of chains = 1 
Sample size per chain = 251 

1. Empirical mean and standard deviation for each variable,
   plus standard error of the mean:

               Mean     SD Naive SE Time-series SE
(Intercept)  1.0736 2.8674  0.18099        0.18099
IDB          -0.9188 0.1922  0.01213        0.01213

2. Quantiles for each variable:

              2.5%    25%     50%     75%   97.5%
(Intercept) -4.952 -0.773  1.1059  3.0165  6.4824
IDB         -1.303 -1.048 -0.9284 -0.7679 -0.5795

the result doesn't like is very stable as it keep changing every time I run it.

Moreover, to find the final offset for the ID B I need to add the (Intercept) Mean to the IDB Mean, how does it work for the SD?

Would it be better to run the spLM formula separately for the two IDs (with y~1 instead of y~ID)?

Thanks

Solution

I am unclear what you mean by "fit an offset model y ~ 1". When I read this, I think you want a model that only has an intercept, but reading further it suggests you want a model where you can estimate the mean for both groups which can be done using

y ~ 0+ID # manually remove the intercept,

To answer your questions:

the result doesn't like is very stable as it keep changing every time I run it.

You are not using very many iterations. Try running with more iterations. With enough iterations the results should be stable.

Moreover, to find the final offset for the ID B I need to add the (Intercept) Mean to the IDB Mean, how does it work for the SD?

Again, I'm not sure what you mean by offset, but if you mean you want the difference in means between group A and group B, this this is exactly what you have in the line beginning with IDB. That is, -0.9188 is the estimated difference in means between group B and group A, i.e. group B's mean is estimated to be 0.9188 smaller than group B's mean, and the SD is the posterior standard deviation.

If you are interested in group B's mean, then you are correct that you must add the (Intercept) to the IDB, but you cannot simply add the SDs. You have two options here: 1) use an appropriate design matrix (X in the above code) that directly obtains your desired parameter estimates or 2) obtain the MCMC samples and calculate the sum of the (Intercept) and IDB parameters for each iteration and then take means and standard deviations of these sums.

Would it be better to run the spLM formula separately for the two IDs (with y~1 instead of y~ID)?

If you ran them separately, then you would be estimating the spatial parameters separately. If the spatial parameters are different in the two different groups, running them separately makes a lot of sense. If they are the same (or similar) then it probably makes more sense to fit the two groups together so you can "borrowing information" about the spatial parameters between these two groups.