I have noticed 3 things:
If s[K] has the initial observations per country, why they are not interacting with y[N] directly?
int s[K];
The second thing is that you are trying to do like an integer auto-regressive model pus new infections so maybe you could use info from this post Log-linear auto-regressive
y[n] ~ poisson(y[n-1] + new_cases);
pos = pos + s[k]; # s[k] not affecting y[n] directly!
And third, the new cases of your model are sampled using the whole information in your process
at last the time, why don’t you use the information of the difference between one time and its last? something like this:
new_cases ~ binomial(y[n] - y[n1], transmission_prob);
Including the recommendations above of not using Poisson directly (you could use log_poisson instead) you should check this post:
In some parts of it, they talk of how to use poisson_glm models and mixed with auto-regressions to simulate the dynamics. And that could be a way of mixing s[k] and y[n-1] at your model. And you should also look for Poisson process your model looks like one.