In measure theory it is possible to extend (for functions from $\mathbb{R}^n$ to $\mathbb{R}$ and integrating wrt Lebesgue measure) one of the direction of the FTC as follows:
For a smooth function of a single variable FTC tells us that
$$f(x)=\frac{d}{dx}\int_a^xf(y)dy=\lim_{r\to 0^+}\frac{1}{2r}\int_{x-r}^{x+r}f(y)dy=\lim_{r\to 0^+}\frac{1}{\
|B_r(x)|}\int_{B_r(x)}f(y)dy$$ where $B_r(x)$ is the ball centered at $x$ and radius $r$ (in this case the ball is an interval) and $|B_r(x)|$ is its measure.
The equation above still holds if $f:\mathbb{R}^n \to \mathbb{R}$ and we use the $\mathbb{R}^n$-Lebesgue measure $d \lambda$ instead of $dy$, provided $f$ is locally integrable (i.e. integrable on every compact set) and that's precisely what Lebesgue differentation theorem states:
$$f(\vec{x})=\lim_{r\to 0^+}\frac{1}{\
|B_r(\vec{x})|}\int_{B_r(\vec{x})}f(\vec{y})d\lambda$$
Now, I was wondering if there is some way to express that using derivatives in $\mathbb{R}^n$ (just as in the one-dimensional case). At first I tried substituting $\frac{d}{dx}$ with the gradient but that doesn't make sense because $f$ is a real-valued function and the gradient is a vector-valued function, any other thoughts?
There is another version of the Lebesgue differentiation theorem, that lets you use e.g. the family of cubes $C_r(\vec x) = [x_1-r,x_1+r]\times\dots\times{}[x_n-r,x_n+r]$ (with $\vec x\in\mathbb R^n$) instead of the balls. Now, since you want to see some kind of derivative you need to assume that $f$ is continuous, not just integrable. You can now consider by keeping Fubinis theorem in mind:
\begin{align*} f(\vec x) &= \frac d{d x_n} \int_0^{x_n} f(x_1,\dots,x_{n-1},y_n)\,d\lambda_1(y_n) \\ &= \frac d{d x_{n-1}}\frac d{d x_n} \int_0^{x_{n-1}}\int_0^{x_n} f(x_1,\dots,x_{n-2},y_{n-1},y_n)\,d\lambda_1(x_n)\,d\lambda_1(x_{n-1}) \\ &= \frac d{d x_1}\dots\frac d{d x_n} \int_0^{x_1}\dots\int_0^{x_n} f(\vec y)\,d\lambda_1(y_n)\dots d\lambda_1(y_1) \\ \Biggl(&= \frac{d^n}{\prod_{k=1}^n d x_k} \int_{\prod_{k=1}^n [0,x_k]} f(\vec y)\,d\lambda_n(\vec y)\Biggr)\\ &= \frac d{d x_1}\dots\frac d{d x_{n-1}} \lim_{r_n\to0+}\frac1{2r_n}\int_0^{x_1}\dots\int_0^{x_{n-1}}\int_{x_n-r}^{x_n+r} f(\vec y)\,d\lambda_n(y_n)\dots d\lambda_1(y_1) \\ &= \lim_{r_1\to0+}\dots\lim_{r_n\to0+}\frac1{2r_1}\dots\frac1{2r_n}\int_{x_1-r}^{x_1+r}\dots\int_{x_n-r}^{x_n+r} f(\vec y)\,d\lambda_1(y_n)\dots d\lambda_1(y_1) \\ &= \lim_{r_1\to0+}\dots\lim_{r_n\to0+}\frac1{2r_1}\dots\frac1{2r_n}\int_{C_r(\vec x)} f(\vec y)\,d\lambda_n(\vec y) \\ &= \lim_{r\to0+}\frac1{(2r)^n}\int_{C_r(\vec x)} f(\vec y)\,d\lambda_n(\vec y) = \lim_{r\to0+}\frac1{|C_r(\vec x)|}\int_{C_r(\vec x)} f(\vec y)\,d\lambda_n(\vec y) \\ \end{align*}
In the second to last line, we see that, since the $n$-fold limit exists (uniformly), the limit with all $r_k$ going at an equal rate $r_k=r\to0$ exists, as well. You can see that
$$f(\vec x)= \lim_{r\to0+}\frac1{|C_r(\vec x)|}\int_{C_r(\vec x)} f(\vec y)\,d\lambda_n(\vec y)= \frac{d^n}{\prod_{k=1}^n d x_k} \int_{\prod_{k=1}^n [0,x_k]} f(\vec y)\,d\lambda_n(\vec y)$$
involves a derivative of order $n$ rather than a gradient.