Work in an isothermal process is expressed as $$W=-\int_{V_0}^{V}\frac{nRT}{V}dV=-nRT\int_{V_0}^{V}\frac{dV}{V}=nRT\ln\left(\frac{V_0}{V}\right)$$
while work in an adiabatic process is expressed as $$W=\frac{P_0V_0}{\gamma - 1}\left[\left(\frac{V_0}{V}\right)^{\gamma -1}-1\right]=\frac{nRT_0}{\gamma - 1}\left[\left(\frac{V_0}{V}\right)^{\gamma -1}-1\right]$$
Is it possible to mathematically prove that expansion work done in an adibatic process is less than that of an isothermal process?
As Poutnik states in the comments, since the integrand is the same, if we can prove that $p_\mathrm{isothermal}(V) > p_\mathrm{adiabatic}(V) \; \forall V$, that would be enough. Since pressure is always a positive physical magnitude, we don't run into problems with justifying cases when the integrand has a negative value.
I will only speak about molar volumes, i.e. $V$ denotes in fact $V/n$, so that the equations are more concise.
Isothermal evolution
In the case of $T = T_0$, the integrand is given directly by the ideal gas law. We modify this a bit for convenience to be seen later, by just dividing each side by the initial pressure $p_0$ \begin{align} \require{cancel} p &= \frac{RT}{V} \\ \frac{p}{p_0} &= \left(\frac{1}{p_0}\right)\frac{RT}{V} \\ \frac{p}{p_0} &= \left(\frac{V_0}{\cancel{RT_0}}\right)\frac{\cancel{RT}}{V} \to \boxed{p = \left(\frac{V_0}{V}\right)p_0} \tag{1} \end{align}
Adiabatic evolution
We can obtain the integrand by manipulation of the 1st law \begin{align} \mathrm{d}U &= \mathrm{d}Q + \mathrm{d}W \tag{$\mathrm{d}Q = 0$} \\ \mathrm{d}U &= \mathrm{d}W \\ C_V\mathrm{d}T &= -p\mathrm{d}V \\ C_V\mathrm{d}T &= -\frac{RT}{V}\mathrm{d}V \\ \int_{T_0}^T \frac{C_V}{RT'} \mathrm{d}T' &= -\int_{V_0}^V\frac{\mathrm{d}V'}{V'} \tag{2} \end{align} We can arrive at the adiabatic function only if we consider that $C_V$ is almost independent of the temperature in the range of integration. Otherwise, the analysis becomes more complicated. We continue with Eq. (2) \begin{align} \frac{C_V}{R}\int_{T_0}^T \frac{\mathrm{d}T'}{T'} &= -\int_{V_0}^V\frac{\mathrm{d}V'}{V'} \\ \frac{C_V}{R} \ln\frac{T}{T_0} &= -\ln\frac{V}{V_0} \qquad\qquad {\text{(Apply ideal gas law twice in $T/T_0$})} \\ \frac{C_V}{R} \ln\frac{pV}{p_0V_0} &= -\ln\left(\frac{V}{V_0}\right) \\ \frac{C_V}{R} \ln\frac{p}{p_0} + \frac{C_V}{R} \ln\frac{V}{V_0} &= -\ln\left(\frac{V}{V_0}\right) \\ \ln\frac{p}{p_0} &= -\left(1 + \frac{R}{C_V}\right)\ln\left(\frac{V}{V_0}\right) \qquad\qquad (\text{$C_V + R = C_p$}) \\ \ln\frac{p}{p_0} &= -\frac{C_p}{C_V}\ln\left(\frac{V}{V_0}\right) \\ \ln\frac{p}{p_0} &= \ln\left(\frac{V}{V_0}\right)^{-C_p/C_V} \\ \frac{p}{p_0} &= \left(\frac{V}{V_0}\right)^{-C_p/C_V} \to \boxed{p = \left(\frac{V_0}{V}\right)^{C_p/C_V}p_0} \tag{3} \end{align}
Now we put Eq. (1) and Eq. (3) in front \begin{align} \boxed{p = \left(\frac{V_0}{V}\right)p_0} && \boxed{p = \left(\frac{V_0}{V}\right)^{C_p/C_V}p_0} \tag{4,5} \end{align} Note that $V_0$ and $p_0$ are constants that denote the starting point of both evolutions.
Since $C_p/C_V > 1$, for every molar volume $V$ the corresponding pressure $p$ will be lower in the adiabatic case. This means that for every molar volume, the adiabatic curve lies below the isothermal curve. Therefore, we have proven that the work done in the adiabatic evolution is lower than the isothermal evolution.
Beneath we have some curves in dimensionless coordinates where the adiabatic evolutions are always below the isothermal evolution. We abuse a bit the notation because $C_p/C_V = 1$ is physically incorrect, but mathematically Eqs. (4) and (5) coincide.
Let:
$W_\mathrm{T}=\mathrm{isothermal\;work}$
$W_\mathrm{A}=\mathrm{adiabatic\;work}$
$y=\frac{W_\mathrm{T}}{W_\mathrm{A}}$
$x=\frac{V_\mathrm{o}}{V}$
$\lambda=\gamma-1$
Dividing both work expressions:
$$y=\frac{W_\mathrm{T}}{W_\mathrm{A}}=\frac{\require{cancel}\cancel{nRT_\mathrm{o}}\ln{\frac{V_\mathrm{o}}{V}}}{\frac{\cancel{nRT_\mathrm{o}}}{\gamma-1}\left[\left(\frac{V_\mathrm{o}}{V}\right)^{\gamma-1}-1\right]}=\frac{(\gamma-1)\ln{\frac{V_\mathrm{o}}{V}}}{\left(\frac{V_\mathrm{o}}{V}\right)^{\gamma-1}-1}=\frac{\ln{\left(\frac{V_\mathrm{o}}{V}\right)}^{\gamma-1}}{\left(\frac{V_\mathrm{o}}{V}\right)^{\gamma-1}-1}=\frac{\ln{x}^{\lambda}}{x^{\lambda}-1}$$
Since we're dealing with an expansion:
$$0<V_\mathrm{o}<V\;\;\implies\;\;0<x<1$$
Since $C_\mathrm{P}$ is always greater than $C_\mathrm{V}$, and the highest possible value of $\gamma$ is approximately 1.66:
$$1<\;\;\gamma=\frac{C_\mathrm{P}}{C_\mathrm{V}}\leq\;\;1.66\;\;\implies\;\;0<\lambda\leq0.66$$
For example, if we plot $y$ vs $x$ with $\lambda=0.25\;$, $\lambda=0.50\;$, and $\lambda=0.66\;$, while respecting the domain of $x$:
We can observe that, regardless of the value of $\lambda$ within its allowed values, $y$ is always greater than 1, which means the isothermal expansion work is always greater than the adiabatic expansion work:
$$y>1\;\;\implies\;\;\frac{W_\mathrm{T}}{W_\mathrm{A}}>1\;\;\implies\;\;W_\mathrm{T}>W_\mathrm{A}$$
