The trick here is to notice that the Witten index for a finite temperature $\beta$ is given by
$$\text{Tr}\left\{(-1)^Fe^{-\beta H}\right\}=\int_{\text{PBC}}\mathcal{D}\phi\mathcal{D}\overline\psi\mathcal{D}\psi\,e^{-S},$$
where the boundary conditions are on a circle of circumference $\beta$. Next, we know that the Witten index is independent of the temperature (it computes the Euler characteristic of the Riemannian manifold), and so we can take the $\beta\to 0$ limit of this expression. In this case, all non-constant modes of the fields $\phi$ and $\psi$ have energy proportional to $1/\beta$, and thus will be exponentially suppressed in the $\beta\to 0$ limit. Thus, the path integral will localize only to those modes which are constant in time, namely
$$\text{Tr}\,(-1)^F\propto\int_{\mathcal{M}}\mathrm{d}\phi\,\sqrt{g}\int\mathrm{d}\overline{\psi}\,\mathrm{d}\psi\,\exp\left(-\frac{\beta}{2}R_{IJKL}\psi^I\overline{\psi}^J\psi^K\overline{\psi}^L\right),$$
where we have traded out our path integral for a standard integral over constant modes only, and the $\sqrt{g}$ term comes from the integral over non-constant modes in the Guassian limit (a factor of $1/\sqrt{g}$ from the bosonic fields and a factor of $g$ from the fermionic ones). The constant of proportionality can be worked out by being careful with the suppression of non-constant modes and working explicitly with the path integral measure over fourier components. However, this is quite technical.
Now, as a warm-up, if the manifold is $2$-dimensional, we have
$$\text{Tr}\,(-1)^F\propto\int\mathrm{d}^2\phi\,\sqrt{g}\int\mathrm{d}^2\overline{\psi}\mathrm{d}^2\psi\exp\left(-\frac{\beta}{2}R_{IJKL}\psi^I\overline{\psi}^J\psi^K\overline{\psi}^L\right).$$
I will leave it to you to show that, when you bring the Grassmann coordinates down from the exponential and integrate over them, the result is
$$\text{Tr}\,(-1)^F\propto\int\mathrm{d}^2\phi\,\sqrt{g}\,R,$$
where $R$ is the Ricci scalar. Since $\text{Tr}\,(-1)^F=\chi(M)$ is the Euler characteristic, this is exactly the statement of the Gauss-Bonnet theorem, up to a constant of proportionality ($1/4\pi$). The technique for higher dimensions can be worked out in a similar fashion.