I feel guilty for bothering other folks with my problems, but after a few days of struggle, I'm running low on options. I seem to be missing the point on how to use bayesmh for custom models.

I'm attempting a simple example of fitting a Monod growth curve to seven observations.
The form of the model is y = θ₁ x / (θ₂ + x)
My intention is to use mean squared error as the negative log likelihood.

Code:
input x y
28 0.053
55 0.060
83 0.112
110 0.105
138 0.099
225 0.122
375 0.125
end

program monod
    args llvar theta1 theta2
    tempvar sq_err
    generate double `sq_err' = ($MH_y - (`theta1' * $MH_extravars / (`theta2' + $MH_extravars ))) ^ 2
    summarize `sq_err', meanonly
    if r(N) < $MH_n {
        scalar `llvar' = .
        exit
    }
    scalar `llvar' = -1 * r(mean)
end

bayesmh y, ///
    llevaluator(monod, parameters({theta1=0.15} {theta2=50}) extravars(x)) ///
    prior({theta1}, flat) prior({theta2}, flat) ///
    noconstant
The resulting parameter estimates are rather poor. That suggests that there is something fundamentally wrong with my approach. I would be grateful if some kind soul can point out my errors.

Thanks in advance,

Stephen