Interesting, I missed that detail in your original message. However, including the initial score as a predictor does imply the assumption that it is measured without error, which is not true (and you don’t assume so for the second score). Theoretically, this may lead to regression towards the mean, which is why I favor a random intercept model. Or just a linear regression of the pair-wise differences in scores. Then you don’t need to assume zero measurement error for the first score.
The approach you propose for adding a random intercept to you Stan model sounds correct. I simpler alternative is to use the package rstanarm (pre-compiled models) or brms (custom models that need to be compiled). Both packages allow the use of standard R formulae as with lm()