Ok, that definitely helps. So if this is similarity (y) by distance (xx), where distance runs from 0 to something large and similarity is bound (by construction?) to the interval [0,1]. I guess my intuition would be to flip this entire thing upside-down and have y be a decreasing function of xx:
μ = pm.Deterministic('μ', 1 / (1 + (xx/κ)) )
This is just a re-specfication of your model and one could prefer one or the other based on intuition. But it might be useful because it highlight at least a couple of things that weren’t necessarily obvious (to me) in the original specification. I would call this a decay function. So first, κ can be thought of as a decay rate. Second, thinking about it this way seems to suggest parametric and functional alternatives to what you have here. For example, this seems to suggest functional forms such as exponential (y=e^{-xx/κ}). If you wanted to keep the current form, there are modification that seem worth considering. For example, it’s not clear to me whether it is uncontroversial to assume that the expected similarity at distance zero must be zero. Your data doesn’t obviously seem to suggest such a requirement (again, eye-ball test), but I also don’t know how your similarity measure is constructed. If you wanted to estimate this as a parameter, you could do this kind of thing:
y_0 = pm.Beta('y_0', alpha=1, beta=5)
μ = pm.Deterministic('μ', y_0 / (1 + (xx/κ)) )
You could do something similar to estimate the asymptotic similarity expected as xx \rightarrow \infty (though your data seems to suggest this is unnecessary).
The other thing that you may want to try is wrapping your mean in a truncated normal. You had asked about asymmetric distributions, but if the skew is due to “artificial” bounds on the outcome variable, truncation may be a more natural solution.