Is the pz orbital wave function three- or four-dimensional?
The wavefunction covers 3-dimensional space. Formally, it maps each point in space to a complex number. That is, it is a function that takes a point in 3-dimensional space as input and returns a complex number as output: $\Psi:\mathbb{R}^{3}\rightarrow \mathbb{C}$.
There are many ways to represent 3-dimensional space: cartesian coordinates $(x,y,z)$, spherical coordinates $(r,\theta,\phi)$, cylindrical coordinates $(z, r, \theta)$, to name the most common. We have these different ways because functions might look &quot;nice&quot; in one coordinate system and not nice in another.
As it turns out, the wave function for the $p_{z}$ orbital looks &quot;nice&quot; in spherical coordinates, that is, $\Psi(r, \theta, \phi)$. It looks less nice in cartesian coordinates, but you can do the replacement for $\theta$ and $\phi$ as noted in the question. However, people do not usually mix coordinates systems together because they are redundant. For example, no one writes $r\sqrt{x^{2}+y^{2}}$. They write either $r^{2}$ or $x^{2}+y^{2}$.
To conclude, the $p_{z}$ orbital you have written is a function of 3 variables: $r, \theta, \phi$. It just happens not to depend on $\phi$, so it looks like just a function of $r$ and $\theta$. Forget this concept of dimensionality because it's not doing you any favors.