Define a Specific Probability Distribution Object in Matlab

I'm trying to create a specific distribution in Matlab to sample from. I want to model some processes that are distributed as ARG (Autoregressive Gamma).

A process v is said to be ARG(a, b, c) if it is of the form v ~ Gamma(a + zt, c)

where zt, ~ Poisson (bvt-1 /c)

The parameters a, b, c correspond, respectively, to a degree of freedom, a measure of serial dependence, and a scale parameter.

I know MATLAB has Poisson and Gamma distributions built in, but I'm not sure how to create a custom distribution function incorporating both of these.

I want to create a function for sampling from processes that follow ARG dynamics. For example, MATLAb has R = mvnrnd(MU,SIGMA)

to generate random numbers from the normal distribution. How can I create something similar for ARG?