Difference between averaging and ignoring the partial dependencies? This question sparks from model interpretation/visualization. To graph the dependency of a function with >2 arguments, one often needs to ignore or average out some arguments.
Problem set
Hastie, Tibshirani & Friedman (Elements of Statistical Learning, 2001) p.370 define a partial dependence of a mapping $f(X)=f(X_S, X_C)$ to a subset of its arguments $X_S$ by averaging the $X_C$ out:
$$\bar{f}_S (X_S)= \mathrm{E}_{X_C} \, f(X_S, X_C)$$ 
and state explicitly that this function does not describe the effect of $X_S$ on f ignoring the effect of $X_C$ on $f$, but $\tilde{f}_S$ does:
$$ \tilde{f}_S (X_S) = \mathrm{E}( f(X_S, X_C) | X_S) $$
Questions
Am I correct in my understanding:


*

*The first equation starts with an "existing" function $f(X_S, X_C)$ and
then removes the dependency on $X_C$ by averaging over it.

*The second equation describes the concept of creating an approximating $\tilde{f}_S$, such that for each "value" of $X_S$ the approximation $\tilde{f}_S$ is close to $f(X_S, X_C)$.


What is the mathematical difference between the two equations? The first expectation integrates only over the $X_C$ dimensions, whereas the second one integrates over all dimensions, right? 
Is it possible to get an analytical sense of the difference of both? 
Is it possible to calculate the second equation in practice for some dataset?
 A: Analytically is a great question. Here's my go at the problem, maybe between the two of us and others we can solve it correctly. 
Let's think of two types of yhats. Let's call them:
\begin{align} y_{hat}= f(Xs=x_s,Xc=x_c)   \end{align}
And, 
\begin{align} y*= f(Xs,Xc=x_c)   \end{align}
A fundamental difference is, one uses the training set exactly (yhat), while y* manipulates the training set. But analytically, they should be different. Here's my first stab at the problem. 
Partial dependence plots use:
\begin{align}  E_{Xc}(Xs=x_s,Xc=x_c)   \end{align}
\begin{align}  f(Xs_j)_{partial dependence} = \sum_{i=1}^n f(Xs=x_{sj},Xc=x_{ci})/N   \end{align}
This is done for each Xs=xs. Note, it's a very real (and probably common) possibility that values conceived were never shown in the training set with the partial dependence plots. Hence you manipulate the training set in some way. 
Whereas the filtering logic is restricted to using values in the training set. For example,
\begin{align}  E(f(Xs=x_s,Xc=x_c)|X_s=x_s)   \end{align}
\begin{align}  f(Xs_j)_{conditionout}=\sum_{j=1}^ny_{hat}/N_j\ ,where \ X_s=x_{sj}  \end{align}
This may not be correct - yet. 
A: 
What is the mathematical difference between the two equations?

These two equations are completely different, but to understand why to let's make a small journey.
By "Law of the unconscious statistician" you can show that $$E[f(x,y)]=\int_{-\infty} ^{+\infty} f(x,y)\cdot p(x,y) dxdy$$
If you don't love math just skip the following two paragraph

Derivation for p.d.f. of $z=f(w)$, where $f$ is not decreasin function.
First easy derivation for the following. If $z=f(w)$ where $f$ deterministic and $w$ is r.v. If $f(w)$ is non decreasing in $w$, then c.d.f. for Z will have the following form:
$$F_z(z)=P(Z<z)=P(f(w)<z)=P(w<f^{-1}(z))=F_w(f^{-1}(z))$$
And p.d.f. for z will have the following form:
$$f_z(z)=F_z'(z)=f_w(w)|_{w=f^{-1}(z)} \cdot (f^{-1}(z))'_{z}$$

Derivation of LOTUS principle for custom functions for monotonically non decreasing function f()
The equation about $E[f(x,y)]$ is not a definition of Expectation by the way. and it's a theorem that people from STATs called LOTUS. Prove for a simple case:
If  $z=f(w),f(w)'>0, \forall w$ then $p_z(z)=p_w(\phi(z)) \cdot |\phi'(z)|=p_x(\phi(z)) \cdot \phi'(z)$. Where $\phi(z)=f^{-1} (z)$.
It's possible to show that from the definition of c.d.f.
and so: $$E[z]=\int_{-\infty} ^{+\infty} z p_z(z) dz=\int_{-\infty} ^{+\infty} z  p_x(\phi(z)) \cdot \phi'(z) dz=\int_{z=-\infty} ^{z=+\infty} z  p_x(\phi(z)) \cdot d\phi(z)=|z=f(w)|=\int_{w=-\infty} ^{w=+\infty} f(w)  p_x(\phi(f(w))) \cdot d\phi(f(w))=\int_{w=-\infty} ^{w=+\infty} f(w)  p_x(w) \cdot dw$$
Maybe it's possible to show via using identity $f^{-1}(y=f(x))=x\implies \frac{df^{-1}}{dy}\cdot \frac{df}{dx}=1$, bit I didn't use it explicitly.
It's possible to prove that theorem for:

*

*Case when the function is monotonically decreasing. If carry math then the equation will be the same

*If a function has a finite number of increasing and decreasing then the formula for $p_z(z)$ is not correct, it will be more complicated.

*And how to prove formula if a number of intervals where function has an infinite number of increasing/decreasing intervals it maybe even trickier.


Lotus allows us to consider the expectation of function as the expectation of a random variable (r.v.) by itself. But this is one of the definitions of expectation through Reiman integral under the assumption that $\int_{-\infty} ^{+\infty} |f(x,y)|\cdot p(x,y) dx < \infty$, and there are distribution that does not satisfy that.
Now this is a conditional expectation (definition combined with LOTUS):
$$g(y_{fixed})=E[f(x,y)|y=y_{fixed}]=\int_{-\infty} ^{+\infty} f(x,y=y_{fixed})\cdot \color{red}{p(x|y=y_{fixed})} dx$$
Partial (Average) Dependency is defined as:
$$F(y_{fixed})=\int_{-\infty} ^{+\infty} f(x,y_{fixed})\cdot \color{red}{p(x)}dx$$

What is the mathematical difference between the two equations?

So this quantity is defined completely differently.

Is it possible to get an analytical sense of the difference of both?
The analytical equations are demonstrated above. Intuition is that conditional expectation is some quantity that takes into account functional dependency between X and Y. But Partial dependence does not at all consider that.

It's strange - you know the value of "Y", but you don't consider this information. Of course, this feature of that tool - partial dependence is slightly another quantity than the conditional expectation of a function of r.v.

The first expectation integrates only over the XC dimensions, whereas the second one integrates overall dimensions, right?

No. the Domain of integration is the same.

Is it possible to calculate the second equation in practice for some datasets?

In practice, if you know p.d.f. and function you can calculate, if you know these functions lie in some span you can calculate it as well. If you don' know p.d.f and you have only samples - you can use Mote-Carlo sampling for evaluating this integral, if your functions $f$, $p$ are from linear space where you evaluate high. dim integrals precisely you can do precisely.
