Python中的统计模拟

import numpy as np

# Define the function to be evaluated
def function(x):
   return x**2

# Generate random samples from a uniform distribution between 0 and 1
samples = np.random.uniform(0, 1, size=10000)

# Evaluate the function at each sample
values = function(samples)

# Compute the average of function values
mean_value = np.mean(values)

print("Mean value of the function: ", mean_value)

import numpy as np

# original sample
data = [1, 2, 3, 4, 5]

# number of bootstrap samples
n_samples = 1000

# array to store the bootstrap samples standard deviation
std_devs = np.empty(n_samples)

# generate bootstrap samples
for i in range(n_samples):
   sample = np.random.choice(data, size=len(data), replace=True)
   std_devs[i] = np.std(sample)

# calculating the lower as well as upper bound of the confidence interval
alpha = 0.05
lower = np.percentile(std_devs, alpha/2*100)
upper = np.percentile(std_devs, (1-alpha/2)*100)

print(f'Confidence interval: [{lower}, {upper}]')

import numpy as np

# target distribution
mean = 0
std = 1
target = lambda x: 1/(std * np.sqrt(2 * np.pi)) * np.exp(-(x - mean)**2 / (2 * std**2))

# proposal distribution
proposal_mean = 0
proposal_std = 2
proposal = lambda x: 1/(proposal_std * np.sqrt(2 * np.pi)) * np.exp(-(x - proposal_mean)**2 / (2 * proposal_std**2))

# initial state
x = 0

# number of samples
n_samples = 10000

# array to store the samples
samples = [x]

# Metropolis-Hastings algorithm
for i in range(n_samples):
   x_new = np.random.normal(x, proposal_std)
   acceptance_prob = min(1, target(x_new) / target(x))
   if np.random.rand() < acceptance_prob:
      x = x_new
   samples.append(x)

# plot the samples
import matplotlib.pyplot as plt
plt.hist(samples, bins=50, density=True)
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# Parameters
p = 0.5  # probability of coin being a head
T = 10  # defining the number of time steps

# Initial state
x = 1

# Random number generator
np.random.seed(0)

# Array to store the states
states = [x]

# Simple stochastic process simulation
for t in range(T):
    x = 1 if np.random.rand() < p else 0
    states.append(x)

print(states)

# Plot the states
plt.step(range(T+1), states, where='post')
plt.xlabel('Time')
plt.ylabel('State')
plt.ylim(-0.5,1.5)
plt.show()

import numpy as np

# Generate random samples from a uniform distribution between 1 and 6
samples = np.random.randint(1, 7, size=10000)

# Compute the average and the standard deviation of the samples
sample_mean = np.mean(samples)
sample_std = np.std(samples)

print("Sample mean: ", sample_mean)
print("Sample standard deviation: ", sample_std)