I have personally also contemplated this issue, and have come up with a simple solution that is satisfactory, to me at least. I'm sure this can also be found in many textbooks. In general, we have
$$\tag{$\star$} Q=\int \rho\ d\tau$$
because we are considering a three dimensional space. Intuitively, we feel it should be possible to talk about a three dimensional charge distribution in every case. The question is how to conceptualize this when discussing surface, line or point charges. The solution comes in the form of the Dirac delta distribution (or function, depending who you ask).
Let's take a look at an example: Consider a 2-sphere of radius $R$, with some charge distribution $\sigma(\theta,\phi)$ on it. What is the three dimensional charge distribution $\rho(r,\theta,\phi)$ corresponding to this situation? Like I said, we have to use the Dirac delta:
$$\rho(r,\theta,\phi)=\delta(r-R)\sigma(\theta,\phi)$$
Now, $(\star)$ gives us:
$$ Q=\int\rho\ d\tau=\int_{0}^{2\pi}\int_{0}^\pi\int_0^\infty \rho\ r^2\sin\theta\ dr\ d\theta\ d\phi=R^2\int_{0}^{2\pi}\int_{0}^\pi\sigma\ \sin\theta\ d\theta\ d\phi$$
Similarly, when considering a line or point charge, one uses two or three Dirac delta's to describe the distribution in 3-space.
This equation will always give you a volume charge density. One way to see this is that surface charge density and volume charge density have different units - $\mathrm{C/m^2}$ and $\mathrm{C/m^3}$ respectively - and in order for the units to be consistent, $\rho$ has to be the latter. The fact that the equation is written with $\rho$ is a helpful reminder that it is a volume charge density.
Of course, keep in mind that the potential is not $kx^{4/3}$ everywhere. That function only describes the potential within a certain region. You also have to think about what's happening outside that region, and on the boundaries of the region.
If you try solving Poisson's equation $\nabla^2\varphi = -\rho/\epsilon_0$ in region where the potential is not so nicely behaved (as you have to do here, if you think about the boundaries), you might get a solution that involves a delta function. Just to pull an example out of thin air, something like
$$\rho(x, y, z) = \delta(x - L) e^{-y^2 - z^2}$$
That is the signature of a surface charge density being expressed as a volume charge density. $\sigma$ is the part other than the delta function; in general:
$$\rho(x, y, z) = \delta(x - a)\sigma(y, z)$$
so in this purely hypothetical example you could discern that $\sigma(y, z) = e^{-y^2 - z^2}$.
This is consistent with the statement that surface charge densities correspond to discontinuities in the electric field, because remember you can write Poisson's equation as
$$\vec\nabla\cdot\vec{E} = \frac{\rho}{\epsilon_0}$$
When $\vec{E}$ is discontinuous, its derivative is "infinite", and therefore $\rho$ needs to be represented as a product involving a delta function.
Best Answer
Sure, though if you want to be super strict, then $\rho$ shouldn't be written inside an integral. You should simply say $\langle\rho,\varphi\rangle:=\int_S\sigma\varphi\,da$, i.e the value of the functional $\rho$ on the test function $\varphi$ is defined to be the integral on the right. This formulation works in any number of dimensions, not just the case of a 2D surface $S$ embedded in $\Bbb{R}^3$.
For instance, let $M$ be a $k$-dimensional embedded submanifold of $\Bbb{R}^n$ ($0\leq k\leq n$), $d\mu_k$ be induced $k$-dimensional volume measure on $M$, and let $\sigma\in L^1(M;d\mu_k)$. We can define a distribution $\delta_{M,\sigma}$ on $\Bbb{R}^n$ by setting for each test function $\varphi$ on $\Bbb{R}^n$ (i.e $\varphi\in C^{\infty}_c(\Bbb{R}^n)$ is smooth compactly supported), \begin{align} \langle\delta_{M,\sigma}, \varphi\rangle&:=\int_M\sigma\varphi\,d\mu_k. \end{align} We refer to $\delta_{M,\sigma}$ as the single-layer distribution on the submanifold $M$ with density $\sigma$. The higher derivatives of this distribution are known as multiplet-layers (the "layer" terminology comes from classical PDE solving methods of single and double layer potentials for Laplace's/Poisson's equation).
See analysis/PDE texts for more information about this sort of stuff, for example, Dieudonne's Treatise on Analysis Vol III (Chapter 17.10), or Hormander's The Analyis of Linear Partial Differential Operators I (Chapter 6.1), or the classic text by Gelfand and Shilov Generalized Functions Volume I (Chapter 3.1 is all about distributions on $\Bbb{R}^n$ concentrated lower dimensional manifolds).