alstat

case <- function(n = 10, mu = 3, sigma = sqrt(5), p = 0.025, rep = 100){
  m <- matrix(NA, nrow = rep, ncol = 4)
  for(i in 1:rep){
    norm <- rnorm(mean = mu, sd = sigma, n = n)
    xbar <- mean(norm)
    low <- xbar - qnorm(p = 1 - p) * (sigma/sqrt(n))
    up <- xbar + qnorm(p = 1 - p) * (sigma/sqrt(n))
  
    if((mu > low) && (mu < up)){
      rem <- 1

alstat

$Matrix
            [,1]      [,2]     [,3] [,4]
  [1,] 3.8195089 2.4336051 5.205413    1
  [2,] 4.6144773 3.2285735 6.000381    0
  [3,] 2.0731334 0.6872296 3.459037    1
  [4,] 3.6641014 2.2781976 5.050005    1
  [5,] 3.8393909 2.4534871 5.225295    1
   :       :         :        :        :
   :       :         :        :        :
 [95,] 1.9025146 0.5166108 3.288418    1

alstat

import numpy as np
import scipy.stats as ss

def case(n = 10, mu = 3, sigma = np.sqrt(5), p = 0.025, rep = 100):
    m = np.zeros((rep, 4))
    
    for i in range(rep):
        norm = np.random.normal(loc = mu, scale = sigma, size = n)
        xbar = np.mean(norm)
        low = xbar - ss.norm.ppf(q = 1 - p) * (sigma / np.sqrt(n))

alstat

system.time(case())

#Output
   user system elapsed
  0.008  0.000   0.008

alstat

import time
t0 = time.time()
case()
time.time() - t0

#Output
0.08859896659851074

alstat

system.time(case(rep = 100000))

#Output
   user system elapsed
  7.076  0.000   7.066

alstat

import time
t0 = time.time()
case(rep = 100000)
time.time() - t0

#Output
64.88077402114868

alstat

import numpy as np
import scipy.stats as ss

def case2(n = 10, mu = 3, sigma = np.sqrt(5), p = 0.025, rep = 100):
    scaled_crit = ss.norm.ppf(q = 1 - p) * (sigma / np.sqrt(n))
    norm = np.random.normal(loc = mu, scale = sigma, size = (rep, n))
    
    xbar = norm.mean(1)
    low = xbar - scaled_crit
    up = xbar + scaled_crit

alstat

case2 <- function(n = 10, mu = 3, sigma = sqrt(5), p = 0.025, rep = 100){
  scaledCrit <- qnorm(p = 1 - p) * (sigma/sqrt(n))
  norm <- matrix(data = rnorm(mean = mu, sd = sigma, n = n*rep), ncol = n, nrow = rep)

  xbar <- rowMeans(norm)
  low <- xbar - scaledCrit
  up <- xbar + scaledCrit

  rem <- (mu > low) & (mu < up)
  m <- cbind(xbar, low, up, rem)

alstat

case3 <- function(n = 10, mu = 3, sigma = sqrt(5), p = 0.025, rep = 100){
  xbar <- rowMeans(
    matrix(rnorm(mean = mu, sd = sigma, n = n*rep),nrow=rep)
  )

  q <- qnorm(p = 1 - p) * (sigma/sqrt(n))
  m <- data.frame(
    xbar,
    low=xbar+q,
    high=xbar-q,