The average of \(mn\) values \(f(x_{ij}^*, y_{ij}^*)\) of a function of two variables can be expressed as a double sum:
\begin{equation*}
\text{Avg} _{mn} = \frac{1}{mn} \sum_{j=1}^n \sum_{i=1}^m f(x_{ij}^*, y_{ij}^*)
\end{equation*}
If we take the limit as \(m\) and \(n\) go to infinity, we obtain what we can define as the average value of \(f\) over the region \(R\)
Although the double sum expression for \(\text{Avg}_{mn}\) appears to be missing an essential component to recognize it as a double Riemann sum that corresponds to a double integral, we can do some algebra to adjust for that. First, note that
\begin{equation*}
\Delta x = \frac{b-a}{m} \quad \quad \text{ and } \quad \quad \Delta y = \frac{d-c}{n}\text{.}
\end{equation*}
Thus,
\begin{equation*}
\frac{1}{mn} = \frac{\Delta x \cdot \Delta y}{(b-a)(d-c)} = \frac{\Delta A}{\text{area of }R},\text{.}
\end{equation*}
Then, the average value of the function \(f\) over \(R\text{,}\) which we denote by \(f_{\operatorname{AVG}(R)}\text{,}\) is given by
\begin{align*}
f_{\operatorname{AVG}(R)} \amp = \lim_{m,n \to \infty} \frac{1}{mn} \sum_{j=1}^n \sum_{i=1}^m f(x_{ij}^*, y_{ij}^*)\\
\amp = \lim_{m,n \to \infty} \frac{1}{A(R)} \sum_{j=1}^n \sum_{i=1}^m f(x_{ij}^*, y_{ij}^*) \cdot \Delta A\\
\amp = \frac{1}{\text{area of }R} \iint_R f(x,y) \, dA.
\end{align*}
Therefore, the double integral of \(f\) over \(R\) divided by the area of \(R\) is the average value of the function \(f\) on \(R\text{.}\)
As an additional note, if
\(f(x, y) \geq 0\) on
\(R\text{,}\) then we can interpret this average value of
\(f\) on
\(R\) as the height of the box with base
\(R\) that has the same volume as the volume of the surface
\(z=f(x,y)\) over
\(R\text{.}\)