# Expression replacement with variable exponents

I have an expression that consists of a large sum of variables of powers like:

X[1]^(2p) X[2]^2 + X[1]^2 X[2]^(2p) + X1[1]^(1+p) X[2]^(1+p) ...


Now I would like to replace X[k]^p with sigma^p (p-1)!!.

While I can use the list index as a variable I cannot use the exponent as variable. For example, I cannot use

expr /. X[k_]^(l_) -> sigma^l (l-1)!!


Instead I have to list all options manually:

expr /. X[k_]^2 X[l_]^(2p) -> sigma^(2+2p) (2p-1)!! /. X[k_]^(p+1) X[l_]^(p+1) -> sigma^(2+2p) p!! p!!


Why is this and how can I implement this?

• I assume that the $\sigma^p (p-1)!!$ comes from if the $X[i]$ have independent normal distributions, then the even central moments are $\sigma^p (p-1)!!$ (with the odd central moments being zero). If you are interested in the expectation of expr, then using Expectation[expr, {X[1] \[Distributed] NormalDistribution[\[Mu], \[Sigma]], X[2] \[Distributed] NormalDistribution[\[Mu], \[Sigma]]}] will avoid the necessity of complex substitution rules. – JimB Jan 30 at 0:16

For the truncated example given, both approaches give the same result.

expr = X[1]^(2 p) X[2]^2 + X[1]^2 X[2]^(2 p) + X[1]^(1 + p) X[2]^(1 + p);

rule = X[k_]^(l_) -> sigma^l (l - 1)!!;

expr /. rule

(* sigma^(2 + 2 p) (p!!)^2 + 2 sigma^(2 + 2 p) (-1 + 2 p)!! *)

expr /. X[k_]^2 X[l_]^(2 p) -> sigma^(2 + 2 p) (2 p - 1)!! /.
X[k_]^(p + 1) X[l_]^(p + 1) -> sigma^(2 + 2 p) p!! p!!

(* sigma^(2 + 2 p) (p!!)^2 + 2 sigma^(2 + 2 p) (-1 + 2 p)!! *)

% == %%

(* True *)


However, neither approach deals with terms of the form X[1] X[2]^(2p) or X[1]^(2p) X[2]

X[1] X[2]^(2 p) + X[1]^(2 p) X[2] /. rule

(* sigma^(2 p) (-1 + 2 p)!! X[1] + sigma^(2 p) (-1 + 2 p)!! X[2] *)


Presumably, you want to handle this with a revised rule that handles a default power. See Optional and Default Arguments

rule2 = X[k_]^(l_.) -> sigma^l (l - 1)!!;

X[1] X[2]^(2 p) + X[1]^(2 p) X[2] /. rule2

(* 2 sigma^(1 + 2 p) (-1 + 2 p)!! *)