Estimating Distribution of Transformed Variables

You’re getting the error because V is a numpy array, but V0 (and all the results of your other operations) are pytensor symbolic tensors. Numpy doesn’t know what to do with those.

Since you’re doing recursive computation, you’re going to need to use a special pytensor operation called scan. You will need to use that instead of a basic python loop. Docs are here, examples of making a scan-based model are here and here.

Since your function is autoregressive, you also get the answer to your distribution question en passant: you will use a pm.CustomDist, with the dist argument as a function that returns the result of your scan, like in the two linked examples. PyMC will automatically infer the logp of your graph and give you the correct distribution.