2

I have a small toy script in Mathematica that I am trying to use to evaluate the pdf of $Y$ where $Y=X^2$, and $X$ is uniformly distributed in $[0,a]$.

The script is,

assum = {a > 0}; 
Subscript[P, X][x_] := If[Inequality[0, Less, x, LessEqual, a], 1/a, 0];

Y[x_] := x^2

Subscript[F, Y][y_] := Integrate[If[Y[x] <= y, 1, 0]*Subscript[P, X][x], 
   {x, -Infinity, Infinity}]

Subscript[P, Y][y_] := D[Subscript[F, Y][y], y]

It shows the correct form of the pdf:

In[192]:=
Simplify[Subscript[P, Y][y], assum]

Out[192]=
Piecewise[{{1/(2*a*Sqrt[y]), y > 0 && a^2 >= y}}, 0]

However, it doesn't evaluate the expression correctly:

In[194]:=
Simplify[Subscript[P, Y][a*(a/2)], assum]

During evaluation of In[194]:= General::ivar:a^2/2 is not a valid variable. >>

Out[194]=
D[1/Sqrt[2], a^2/2]

The session, in pretty printing looks like this:

enter image description here

highBandWidth
  • 403
  • 2
  • 10

2 Answers2

7

$P[y]$ is defined to be $\partial_yF[y]$. If you now call $P[2]$ for example, this is translated to $\partial_2F[2]$, which of course does not make any sense.

One way to get around this behavior is not taking the derivative with respect to y, but with respect to the first argument instead. The following example assigns f[y] to be the derivative of g[y]:

g[y_] := y^2
f[y_] := Derivative[1][g][y]
f[x]

2 x

Instead of Derivative[1][g][y] you could also have used the shorthand notation g'[y]. Read the explicit long version I used as an operator applied to multiple things: Derivative[1] takes a function and calculates its first derivative with respect to the first argument. Derivative[1][g] is that derivative (as a pure function), and Derivative[1][g][y] is that pure function applied to a value y.

The script above has one problem however: Every time you call f[x], it re-calculates the derivative, which can take some time if your function is more complicated and you need a lot of data points. If you use functions instead of patterns (i.e. no :=) you can also get around this problem:

g = #^2 &
f = g'
f[x]

#1^2 &
2 #1 &
2 x
David
  • 14,911
  • 6
  • 51
  • 81
  • Thanks! I didn't realise := wasn't the right way to define functions. – highBandWidth Mar 20 '12 at 05:54
  • @highBandWidth err, both = and := are right way to define functions. There is a big difference between them however, and you should learn the difference. See for example http://mathematica.stackexchange.com/questions/2639/what-is-x-x-trickery and SetDelayed: http://reference.wolfram.com/mathematica/ref/SetDelayed.html – tkott Mar 20 '12 at 10:53
  • 1
    @tkott That's not really true. = and := define pattern-based replacement rules, which may look like functions in many cases, but are handled differently internally. Defining actual functions is done using Function. Related: http://mathematica.stackexchange.com/questions/704/functions-vs-patterns – David Mar 20 '12 at 17:54
  • @David Sorry you're absolutely right. I just used the wrong terminology. – tkott Mar 20 '12 at 18:18
2

Or you can just keep your code the way it is (with the SetDelayed), but force the derivation to be executed at the definition by using Evaluate:

Subscript[P, Y][y_] := Evaluate[D[Subscript[F, Y][y], y]]

This gives you the correct result:

In[107]:= Simplify[Subscript[P, Y][a*(a/2)], assum]
Out[107]= 1/(Sqrt[2] a^2)

This would be the most robust way to do it (at least it is how the Mathematica 8 help does it).

freddieknets
  • 1,085
  • 8
  • 16
  • This is how I would have done it, until I learned that one can also just use Subscript[P, Y][y_] = D[Subscript[F, Y][y], y]. The only difference to the original form is that I don't use SetDelayed but Set, so that the y is used literally and not replaced by its value at the time the "function" is called. This way is just a little shorter than using Evaluate. – Jens Jun 10 '12 at 18:16
  • True, but in general I try to avoid the use of Set, as it can give rise to difficult-to-debug naming conflicts. – freddieknets Jun 10 '12 at 21:41