Edit: I have solved this. I was incorrectly simulating values for the dependent variable and not using a good summary statistic. I have changed my algorithm to:

  1. Simulate y value using simulated B0, B1, B2. Random normal around simulated sigma.
  2. Calculate summary statistic as sum of squared differences between simulated y and y from observed data at the same x

I am practicing ABC using a regression model. x is estimated as $\beta_0$ + $\beta_1$*predictor + $\beta_2$*predictor$^2$.

Using lm() in R I am getting this:

data.lm <- lm(dependent ~ predictor + I(predictor^2), data = data)

lm(formula = dependent ~ predictor + I(predictor^2), data = data)

    Min      1Q  Median      3Q     Max 
-2.8411 -0.9694  0.0017  1.0181  3.3900 

                 Estimate Std. Error      t value Pr(>|t|)    
(Intercept)      -7.5555495    1.4241091  -5.305  1.69e-05 ***
Predictor        1.2716937     0.0757321  16.792  3.99e-15 ***
I(Predictor ^2)  -0.0145014    0.0008719  -16.633 4.97e-15 ***
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.663 on 25 degrees of freedom
Multiple R-squared:  0.9188,    Adjusted R-squared:  0.9123 
F-statistic: 141.5 on 2 and 25 DF,  p-value: 2.338e-14

I am implementing the algorithm as follows:

  1. Set priors for $\beta_0$, $\beta_1$, $\beta_2$ each as $\mu$ = 0, $\sigma$ = 10
  2. Simulate values for each beta using rnorm(mean, sd)
  3. Set a linear model using the simulated y values and x values from the data

    n <- length(data_predictor)
    y <- rnorm(n,b0_r[i] + b1_r[i]*x + b0_r[i]*x^2, sigma_r[i])
    sim.lm <- lm(y~data_predictor + I(data_predictor^2))
  4. Calculate the sum of squared distance from the simulated model's error to the error of the model using the data

  5. If distance[i] < epsilon (currently 5), accept and put each beta into accepted values arrays

The problem with this is I am getting inaccurate simulated models. The predicted y values are either extremely negative or positive, which gives a high error in the simulated linear model. I know that ABC is heavily dependent on the priors, so do my priors for $\beta_0$, $\beta_1$, $\beta_2$ need to be close to the coefficients from the linear model using the data?

Is there a correct way to implement ABC on a simple linear model?


closed as unclear what you're asking by Michael Chernick, kjetil b halvorsen, mdewey, Peter Flom Oct 1 '17 at 12:31

Please clarify your specific problem or add additional details to highlight exactly what you need. As it's currently written, it’s hard to tell exactly what you're asking. See the How to Ask page for help clarifying this question. If this question can be reworded to fit the rules in the help center, please edit the question.

  • $\begingroup$ what is distance[i] ? $\endgroup$ – Xi'an Oct 3 '17 at 8:26

Both your standard deviation and epsilon are incredibly high. I would shrink them down and see how that affects your distributions. Maybe try setting them both to one and see if your distributions look more like what you want? These are both tuning parameters that you will need to tweak until things look good.

Take caution, though, because too small of a standard deviation will not explore the parameter space adequately, whereas too small of an epsilon will result in very few accepted parameters!


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