Rayat Shikshan Sanstha’s
Yashavantrao Chavan Institute of Science, Satara (Autonomous)
Department of Statistics
MSP 406: Planning and Analysis of Industrial Experiments
Monte Carlo Simulation
Practical No: 16 Date:
Name: Bagadi Shashikant Pundlik Roll No: 232
Q1
Determine the probability of type I error (false alarm rate) and the performance of the control
chart
when the process mean shifts from 50 to 55. Use Monte Carlo simulation to evaluate the
performance of an 𝑋-bar control chart for monitoring the mean of a process. Assume that our
process is normally distributed with: Mean (𝜇) = 50, Standard deviation (𝜎) = 10, Sample size
(𝑛)
= 5.
1. Write the steps in Monte Carlo simulation.
2. Create the control chart:
3. Evaluate the control chart and Compute Control Limits
4. Write SAS code
5. Drown your conclusion
Q2
Using Monte Carlo simulation to estimate the probability that the sample mean of a given
population exceeds a certain value. Suppose that a population that follows a normal
distribution
with a mean (𝜇μ) of 50 and a standard deviation (𝜎σ) of 10. We want to estimate the
probability
that the sample mean of a sample of size 30 exceeds 52.
1. Write the steps in Monte Carlo simulation.
2. Write SAS code
3. Drown your conclusion
Q1
Aim: To determine the probability of type I error (false alarm rate) and the performance of the
control chart when the process mean shifts from 50 to 55.
1. Write the steps in Monte Carlo simulation.
2. Create the control chart:
3. Evaluate the control chart and Compute Control Limits
4. Write SAS code
5. Down your conclusion
Given: Mean (μ) = 50, Standard deviation (σ) = 10, Sample size (n) = 5
Solution:
Steps in Monte Carlo simulation.
1. Define the Problem and Objectives
2. Develop the Mathematical Model.
3. Identify and Characterize Input Variables.
4. Generate Random Inputs.
5. Run the Simulation.
6. Repeat the Simulation
7. Analyze the Results.
8. Interpret the Results.
9. Validate the Model:
Code:
%let mu = 50;
%let sigma = 10;
%let n = 5;
%let N = 5000;
data simulation1;
call streaminit(12345);
do i = 1 to &N;
sample_mean = mean(rand("Normal", &mu, &sigma),
rand("Normal", &mu, &sigma),
rand("Normal", &mu, &sigma),
rand("Normal", &mu, &sigma),
rand("Normal", &mu, &sigma));
output;
end;
run;
/* Compute control limits */
proc means data=simulation1 noprint;
var sample_mean;
output out=limits mean=mean_xbar std=std_xbar;
run;
data limits1;
set limits;
ucl = mean_xbar + 3 * std_xbar/sqrt(&n);
lcl = mean_xbar - 3 * std_xbar/sqrt(&n);
run;
proc print data=limits1;
var mean_xbar std_xbar ucl lcl;
run;
/* Evaluate Type I error rate */
data type1_error1;
set simulation1;
if sample_mean > 63.41 or sample_mean < 39.59 then type1_error = 1;
else type1_error = 0;
run;
proc freq data=type1_error1;
tables type1_error;
run;
/* Evaluate the control chart when the mean shifts from 50 to 55 */
data simulation_shift1;
call streaminit(12345); /* Set random seed for reproducibility */
do i = 1 to &N;
sample_mean = mean(rand("Normal", 55, &sigma),
rand("Normal", 55, &sigma),
rand("Normal", 55, &sigma),
rand("Normal", 55, &sigma),
rand("Normal", 55, &sigma));
output;
end;
run;
proc means data=simulation_shift1 noprint;
var sample_mean;
output out=limits mean=mean_xbar std=std_xbar;
run;
data limits2;
set limits;
ucl = mean_xbar + 3 * std_xbar/sqrt(&n);
lcl = mean_xbar - 3 * std_xbar/sqrt(&n);
run;
proc print data=limits1;
var mean_xbar std_xbar ucl lcl;
run;
/* Calculate the rate of detecting the shift */
data detection_rate;
set simulation_shift1;
if sample_mean > 63.41 or sample_mean < 39.59 then detected = 1;
else detected = 0;
run;
proc freq data=detection_rate;
tables detected;
run;
Output:
Obs mean_xbar std_xbar ucl lcl
1 50.0868 4.40914 50.2738 49.8997
Cumulative Cumulative
type1_error Frequency Percent Frequency Percent
0 149 2.98 149 2.98
1 4851 97.02 5000 100.00
Obs mean_xbar std_xbar ucl lcl
1 50.0868 4.40914 50.2738 49.8997
Cumulative Cumulative
detected Frequency Percent Frequency Percent
0 71 1.42 71 1.42
1 4929 98.58 5000 100.00
Code:
%let mu = 50;
%let sigma = 10;
%let n = 5;
%let N = 5000; /* Number of simulations */
/* Generate random samples and calculate sample means */
data simulation;
call streaminit(12345); /* Set random seed for reproducibility */
do i = 1 to &N;
do j = 1 to &n;
x = rand("Normal", &mu, &sigma);
output;
end;
end;
run;
/* Calculate sample means */
proc means data=simulation noprint;
by i;
var x;
output out=sample_means mean=mean_xbar;
run;
/* Create X-bar control chart */
proc sgplot data=sample_means;
series x=i y=mean_xbar / markers;
refline 50 / lineattrs=(color=red);
refline 36.59 / lineattrs=(color=red);
refline 63.41 / lineattrs=(color=red);
xaxis label="Sample Number";
yaxis label="X-bar";
keylegend / location=inside position=topright across=1;
run;
Conclusion
1. The Type I error rate can be determined by the proportion of sample means outside
the control limits when the process mean is 50.
2. The effectiveness of the control chart in detecting a mean shift from 50 to 55 can be
assessed by the proportion of sample means outside the control limits after the shift.
3. The X- bar control chart with the computed control limits is expected to effectively
monitor the process and detect shifts in the mean, provided the assumptions of
normality and consistent standard deviation hold true.
Q2.
Aim: To estimate the probability that the sample mean of a sample of size 30 exceeds 52.
1. Write the steps in Monte Carlo simulation.
2. Write SAS code
3. Drown your conclusion
Given: Mean (μ)=50, Standard deviation (σ)= 10.
Solution:
Steps in Monte Carlo simulation.
1. Define the Problem and Objectives
2. Develop the Mathematical Model.
3. Identify and Characterize Input Variables.
4. Generate Random Inputs.
5. Run the Simulation.
6. Repeat the Simulation
7. Analyze the Results.
8. Interpret the Results.
9. Validate the Model:
Code:
%let mu = 50;
%let sigma = 10;
%let sample_size = 30;
data monte_carlo;
do i = 1 to 10000
array sample[&sample_size];
do j = 1 to &sample_size;
sample[j] = rand("Normal", &mu, &sigma);
end;
sample_mean = mean(of sample[*]);
if sample_mean > 52 then exceeds = 1;
else exceeds = 0;
output;
end;
run;
proc means data=monte_carlo noprint;
var exceeds;
output out=summary mean=probability;
run;
proc print data=summary;
title 'Probability of Sample Mean Exceeding 52';
var probability;
run;
Output:
Probability of Sample Mean Exceeding 52
Obs probability
1 0.1409
Result:
The probability that the sample mean exceeds 52 is 0.1409
