Can I Expand "Nested" Symbolic Functions in Matlab?

I'm trying to model the effect of different filter "building blocks" on a system which is a construct based on these filters. I would like the basic filters to be "modular", i.e. they should be "replaceable", without rewriting the construct which is based upon the basic filters.

For example, I have a system of filters G_0, G_1, which is defined in terms of some basic filters called H_0 and H_1.

I'm trying to do the following:

    syms z
    syms H_0(z) H_1(z)
    G_0(z)=H_0(z^(4))*H_0(z^(2))*H_0(z)
    G_1(z)=H_1(z^(4))*H_0(z^(2))*H_0(z)

This declares the z-domain I'd like to work in, and a construct of two filters G_0,G_1, based on the basic filters H_0,H_1.

Now, I'm trying to evaluate the construct in terms of some basic filters:

    H_1(z) = 1+z^-1
    H_0(z) = 1+0*z^-1

What I would like to get at this point is an expanded polynomial of z. E.g. for the declarations above, I'd like to see that G_0(z)=1, and that G_1(z)=1+z^(-4).

I've tried stuff like "subs(G_0(z))", "formula(G_0(z))", "formula(subs(subs(G_0(z))))", but I keep getting result in terms of H_0 and H_1.

Any advice? Many thanks in advance.

Edit - some clarifications:

In reality, I have 10-20 transfer functions like G_0 and G_1, so I'm trying to avoid re-declaring all of them every time I change the basic blocks H_0 and H_1. The basic blocks H_0 and H_1 would actually be of a much higher degree than they are in the example here.
G_0 and G_1 will not change after being declared, only H_0 and H_1 will.
H_0(z^2) means using z^2 as an argument for H_0(z). So wherever z appears in the declaration of H_0, z^2 should be plugged in
The desired output is a function in terms of z, not H_0 and H_1.
A workable hack is having an m-File containing the declarations of the construct (G_0 and G_1 in this example), which is run every time H_0 and H_1 are redefined. I was wondering if there's a more elegant way of doing it, along the lines of the (non-working) code shown above.

Solution

This seems to work quite nicely, and is very easily extendable. I redefined H_0 to H_1 as an example only.

syms z

H_1(z) = 1+z^-1;
H_0(z) = 1+0*z^-1;

G_0=@(Ha,z) Ha(z^(4))*Ha(z^(2))*Ha(z);
G_1=@(Ha,Hb,z) Hb(z^(4))*Ha(z^(2))*Ha(z);

G_0(H_0,z)
G_1(H_0,H_1,z)

H_0=@(z) H_1(z);

G_0(H_0,z)
G_1(H_0,H_1,z)