This looks similar to what in ecology is often called an N-mixture model - that might help you find relevant implementations and more details. The trick you need is that you have to marginalise (sum up) over all possible values of ystar
(from y to infinity; in practicality, you will need to assign an upper bound above which you believe the probability of ystar is negligible), weighted by the probability of that value of ystar
under the Poisson likelihood.
The likelihood for a single data point looks like the following, with an upper bound of 10000 for demonstration (this will depend on your data):
int ub = 10000;
vector[ub - y + 1] lp;
for(i in 1:(ub - y + 1)) {
lp[i] = binomial_lpmf(y | i, p) + poisson_log_glm_lpmf(i | x, alpha, beta);
}
target += log_sum_exp(lp);
You can find an example implementation of an N-mixture model like this here: Calculating abundance from marginalized N-mixture model
edit: for more information on marginalisation, the Stan users guide is very good: Latent Discrete Parameters, as is @mbjoseph’s post about occupancy models, which has a great beginner’s guide to the mathematics: Maxwell B. Joseph: A step-by-step guide to marginalizing over discrete parameters for ecologists using Stan