Error in calculation of $\pi$ using Monte Carlo method. So lets say we are trying to calculate value of $\pi$ using MonteCarlo method. By picking random points in a square and measuring their distance from the center and if $k$ points lie inside the circle using the ratio $\frac{k}N$ to calculate the value of $\pi$. 
How do I derive the formula for variance of the calculated value of $\pi$?
This is what I have got so far. 


*

*Say the probability that the point lies inside the circle is $p$. 

*Say we have N random points in the square, of which k lie inside the circle. 

*The probability that k points lie inside the circle is $\binom{N}{k}p^k(1-p)^{N-k}$. 

*The error in estimate of $\pi$ is $e=||c\pi-\frac{N}k||$. 

*So the expected value of error $E(e) = \sum_{k=0}^{k=N}\left\|c\pi-\frac{N}k\right\|\binom{N}{k}p^k(1-p)^{N-k}$?
Is the expectation value of error that I have written correct, how do I evaluate it?
 A: For very large $n$, as in a simulation, a reasonable 95% confidence interval
for the estimate $\hat p = k/n$ of the proportion of the points inside the circle is of the form
$$\hat p \pm 1.96\sqrt{\frac{\hat p(1-\hat p)}{n}}.$$
So the approximate margin of simulation error is $1.96\sqrt{\frac{\hat p(1-\hat p)}{n}}.$ 
We have $Var(\hat p) = Var(k/n) = (1/n^2)Var(k) = p(1-p)/n,$ so
 $SD(\hat p) = \sqrt{\frac{p(1-p)}{n}},$ where 
$p$ is the unknown theoretical fraction. We use the observed value $\hat p$
to get $\widehat{SD}(\hat p),$ the estimated standard deviation of $\hat p.$ 

But be careful! You are not specific about 'a square' and 'the circle', but the usual
simulation is to use a square of area $4$ with vertices at $(-1,-1)$ and $(1,1)$
and the inscribed circle of radius $1$ and area $\pi.$ Thus the ratio of
'accepted' points in the circle to the total points in the square estimates
$\pi/4.$ If you multiply the estimate by $4$ to get an estimate of $\pi,$
then the margin of error is also multiplied by $4.$
Let's try the actual simulation experiment in R statistical software:
> set.seed(1066)  # uae this seed to replicate this simulation
> n = 10^6;  x = runif(n, -1, 1);  y = runif(n, -1, 1)
> r = sqrt(x^2 + y^2)
> p = sum(r < 1)/n;  pi.est = 4*p;  pi.est
[1] 3.140892
> abs(pi.est - pi)  # observed simulation error
[1] 0.0007006536
> me.p = 2*sqrt(p*(1-p)/n); me.p
[1] 0.0008213351
> me.pi = 4*me.p;  me.pi  # 95% est. marg of sim err, using sim results
[1] 0.00328534
q = pi/4
4*1.96*sqrt(q*(1-q)/m)    # exact 95% marg of sim err, using pi=3.14.59...
[1] 0.003218679

The observed simulation error (about 0.001) for estimating $\pi$ 
is within the 95% simulation margin of error (about 0.003) for the estimate of $\pi$,
but not within the 95% margin of simulation error (about 0.0001) for estimating
$\pi/4.$
Here is a plot of 10,000 of the million points.

