Recently I notice Stan new reduce-sum function which could parallelly compute the log-likelihood. I am wondering could we use it in Rstan now or it is still onlyavailable in cmdstanr? Also is it applicable in any of the Bayesian modeling in Stan(Since I think every Bayesian model must have log-lieklihood and thus it seems that reduce-sum could be applied in every Stan script)?

It can be used in a lot of them! I think next release of rstan that uses the new compiler will support `reduce_sum`

in rstan. This intro should help you get started

Yes this intro is very helpful. And I am curious about whether we could use reduce_sum in any of the Stan script since every model (eg hierarchical model) will involve log likelihood function?

As long as things are conditionally independent you can use reduce-sum. That’s for a lot, but - for example - not for Gaussian Processes where you need to invert a big fat covariance matrix in a singe operation.

Thanks for the explanation. So may I ask that do we have any examples about how could we use the reduce-sum in hierarchical models?

Look at the brms vignette for threading.

Another example is Simon Mills’ vignette on using reduce_sum with a multi-species occupancy model (a hierarchical model, in this example parameterized with non-centered random effects), available at https://github.com/SimonCMills/reduce_sum_example