How to implement a formal expectation operator over an unknown distribution?

Question

I need to do some simplification of an expression involving averages over a stochastic variable (in order to verify a long analytical calculation). The easiest way to do that, I figured, were if I could implement an operator which would basically be short-hand for the averaging procedure, with all the appropriate properties. Then of course this operator would be present in the final expression, which is fine, and would enable me to compare easily with my own calculations. So assuming I use x for the stochastic variable, I tried defining av using

av[y_ + z_] := av[y] + av[z]  
av[c_ y_] := c av[y] /; FreeQ[c, x]
  
av[c_] := c /; FreeQ[c, x]  

Then when I write

D[av[x y], y]

I get


av[x]

which is fine, but when I write

D[av[Exp[-x y]], y]

I get

-E^(-x y) y  

instead of -y av[Exp[-x y]] as I want, i.e., the av is removed somehow.

I tried using UpValue for teaching Mathematica that it could interchange differentiation and av, but apparently that is not the problem. I might be going about this entirely the wrong way, but I'd be grateful for any input. Note the builtin Expectation function does not accomplish it either - e.g., it doesn't handle the derivatives as a proper average operator would. For example

h[y_] := Expectation[y, x \[Distributed] pp] (*pp unknown density*)

Then

D[h[Exp[-x y]], x]

gives

(Expectation^(0,1))[E^(-x y),x\[Distributed]pp] (Distributed^(1,0))[x,pp]-E^(-x y) y   
(Expectation^(1,0))[E^(-x y),x\[Distributed]pp]

whereas I wanted

-y h[Exp[-y x]]

(i.e moving the derivative inside the averaging h). Sune

Answer 1

I think I got it.

av /: D[av[f___], x_] := av[D[f, x]]
av[y_ + z_] := av[y] + av[z]
av[c_ y_] := c av[y] /; FreeQ[c, x]
av[c_] := c /; FreeQ[c, x]
D[av[x y], x]

(* y *)

D[av[Exp[-x y]], x]
(* -y av[E^(-x y)] *)

I wasn't using UpValues correctly before.

EDIT: Well, turns out there's still a problem:

D[Log[av[Exp[-b x]]], b]
(* -((E^(-b x) x)/av[E^(-b x)]) *)

instead of

-(av[(E^(-b x) x)]/av[E^(-b x)])

What's going on? Sune

Answer 2

Just define exactly what you want:

av[expr_] := Integrate[f[x] expr, {x, -∞, ∞}];

Here f[x] is undefined function.

D[av[Exp[-x y]], y]

D[Log[av[Exp[-b x]]], b]

P.S. There is a possible bug

f[c_] := c /; FreeQ[c, x]

f'
(* 1 & *)

The condition was lost.