Finding the mean and variance of a distribution

Question

I have tried unsuccessfully to use Mathematica to derive both the expected value (mean) and variance of a particular statistical distribution of interest to me.

The distribution is a discretized version of the Gaussian, and the question is: what is the expected value of $\lfloor X \rfloor$, where $X\sim N(\mu, \sigma)$ and $\lfloor \cdot \rfloor$ denotes the floor function (i.e., in Mathematica, Floor[]) ?

I have tried a direct integration approach (with assumptions) allowing for a general mean and standard deviation:

Integrate[
 Floor[x]*PDF[NormalDistribution[\[Mu], \[Sigma]], 
   x], {x, -\[Infinity], \[Infinity]}, 
 Assumptions -> {\[Mu] \[Element] Reals, \[Sigma] \[Element] 
    PositiveReals}]

and also tried restricting the question to one involving an underlying standard normal distribution:

Integrate[
 Floor[x]*PDF[NormalDistribution[0, 1], 
   x], {x, -\[Infinity], \[Infinity]}]

Equivalently, for the 2nd raw moment:

Integrate[
 (Floor[x])^2 * PDF[NormalDistribution[\[Mu], \[Sigma]], 
   x], {x, -\[Infinity], \[Infinity]}, 
 Assumptions -> {\[Mu] \[Element] Reals, \[Sigma] \[Element] 
    PositiveReals}]

The approach fails in both instances, but the problem is solvable. One solution, for the mean, which uses the Poisson summation formula, is:

$$ \mathbb E\lfloor X \rfloor = \mu - \mathbb E \{X\} = \mu - \frac{1}{2} + \frac{1}{\pi} \sum_{k=1}^\infty e^{-2\pi^2\sigma^2 k^2} \frac{\sin(2\pi k \mu)}{k} $$

and, indeed, it strikes me that using the Poisson summation formula is likely to yield a solution for the question about variance also. However, to use the Poisson summation formula, one needs the Fourier transform of the distribution, but with that I have also been unsuccessful.

What is an appropriate way of using Mathematica to obtain a solution both for the mean, and the variance for my distribution.

Further explanation: I've just seen Bob Hanlon's very instructive numerical-methods solution ... but I'd really like to find a way of reproducing the Poisson summation/Fourier-based analytic solution. One very significant reason for preferring a solution involving the Poisson summation formula (and associated Fourier description) is that it clearly reveals the oscillatory nature of the function, and it is the oscillation that particularly interests me. It's true that I have displayed an analytic description of the mean for the Floor[] function, but there are several other distributions of interest, so I'll like to learn a general Mathematica approach.

Answer 1

Use numeric methods.

Clear["Global`*"]

For the mean,

μ[m_?NumericQ, s_?Positive] := NIntegrate[Floor[x]*
     PDF[NormalDistribution[m, s], x],
    {x, -∞, ∞},
    WorkingPrecision -> 30,
    MinRecursion -> 25,
    MaxRecursion -> 100] //
   N // Quiet

For example,

μ @@@ {{-3, 2}, {0, 1}, {3, 2}}

(* {-3.5, -0.5, 2.5} *)

For the variance,

var[m_?NumericQ, s_?Positive] := 
   NIntegrate[(Floor[x] - μ[m, s])^2*
     PDF[NormalDistribution[m, s], x],
    {x, -∞, ∞},
    WorkingPrecision -> 30,
    MinRecursion -> 25,
    MaxRecursion -> 100] //
   N // Quiet

Examples,

var @@@ {{-3, 2}, {0, 1}, {3, 2}}

(* {4.08333, 1.08333, 4.08333} *)

Answer 2

A shorter way is as follows.

N[Mean[TransformedDistribution[Floor[x], x \[Distributed] NormalDistribution[-3, 2]]], 3]

-3.50

N[Variance[TransformedDistribution[Floor[x], x \[Distributed] NormalDistribution[-3, 2]]], 3]

4.08

Addition. One more way is as follows.

NSum[n*Integrate[ PDF[NormalDistribution[-3, 2], t], {t, n, n + 1}], {n, -Infinity, Infinity}]

-3.5

NSum[n^2*Integrate[PDF[NormalDistribution[-3, 2], t], {t, n, n + 1}], {n, -Infinity,  Infinity}] -NSum[n*Integrate[PDF[NormalDistribution[-3, 2], t], {t, n, n + 1}], {n, -Infinity, Infinity}]^2

4.083333333333332

Answer 3

This is another numerical approach that seems to be much quicker than the other numerical approaches. This simply starts out acknowledging that you have a discrete distribution over the integers. We have

$$\text{Pr}(\lfloor X \rfloor = x) =\Phi \left(\frac{x+1-\mu}{\sigma }\right)-\Phi \left(\frac{x-\mu }{\sigma }\right)=\frac{1}{2} \text{erf}\left(\frac{\mu -x}{\sqrt{2} \sigma }\right)-\frac{1}{2}\text{erf}\left(\frac{\mu -x-1}{\sqrt{2} \sigma }\right)$$

A function that returns the original mean and standard deviation for $X$ and the estimated mean and variance for $\lfloor X \rfloor$ follows:

stats[μ_, σ_] := Module[{n1, n2, mean, variance},
  n1 = Floor[μ - 6 σ];
  n2 = Ceiling[μ + 6 σ];
  If[IntegerQ[μ], mean = μ - 1/2,
   mean = Sum[x (-(1/2) Erf[(-1 - x + μ)/(Sqrt[2] σ)] + 
     1/2 Erf[(-x + μ)/(Sqrt[2] σ)]), {x, n1, n2}]];
  variance = Sum[(x - mean)^2 (-(1/2) Erf[(-1 - x + μ)/(Sqrt[2] σ)] + 
     1/2 Erf[(-x + μ)/(Sqrt[2] σ)]), {x, n1, n2}];
  {μ, σ, mean // N, variance // N}
]

For example:

stats[-3, 2] // AbsoluteTiming
(* {0.0010298, {-3, 2, -3.5, 4.08333}} *)
stats[0, 1] // AbsoluteTiming
(* {0.0004044, {0, 1, -0.5, 1.08333}} *)
stats[3, 2] // AbsoluteTiming
(* {0.0009363, {3, 2, 2.5, 4.08333}} *)

Note that for $\sigma\geq 1$ even quicker approximations for the mean and variance are, respectively:

$$\mu-1/2$$ $$\sigma^2+1/12$$