Discussion Problem

Attempt to complete either 5.9, 5.13, or 5.14 in Rizzo’s Statistical Computing in R. Post your attempts here and discuss.

Problem 5.9

The Rayleigh density [156, (18.76)] is \(f(x)=\frac { x }{ { \sigma }^{ 2 } } { e }^{ { -{ x }^{ 2 } }/{ \left( 2{ \sigma }^{ 2 } \right) } }\) for \(x\geq 0,{ \sigma }>0\). Implement a function to generate samples from a Rayleigh(\(\sigma\)) distribution, using antithetic variables. What is the percent reduction in variance of \(\frac { X+{ X }^{ ' } }{ 2 }\) compared with \(\frac { X_{ 1 }+{ X_{ 2 } } }{ 2 }\) for independent \(X_1,X_2\)?

rayleigh_red <- function(scale, n) {
  rayleigh <- antithetic <- numeric(n)
  for (i in 1:n) {
    U <- runif(n)
    V <- 1 - U
    rayleigh = scale * sqrt(-2 * log(U))
    antithetic = scale * sqrt(-2 * log(V))
    var1 <- var(rayleigh)
    var2 <- (var(rayleigh) + var(antithetic) + 2 * cov(rayleigh, antithetic)) / 4
    reduction <- ((var1 - var2) / var1)
    percent <- paste0(formatC(100 * reduction, format = "f", digits = 2), "%")
  }
  return(noquote(percent))
}
scale = 2
n <- 1000
rayleigh_red(scale, n)
## [1] 96.91%

Problem 5.13

Find two importance functions \(f_1\) and \(f_2\) that are supported on \((1, \infty)\) and are close to \(g(x)=\frac { { x }^{ 2 } }{ \sqrt { 2\pi } } { e }^{ { -{ x }^{ 2 } }/{ 2 } }\) for \(x>1\). Which of your two importance functions should produce the smaller variance in estimating \(\int _{ 1 }^{ \infty }{ \frac { { x }^{ 2 } }{ \sqrt { 2\pi } } { e }^{ { -{ x }^{ 2 } }/{ 2 } } } dx\).

\[f_1={ e }^{ -x }, \quad x \in (1, \infty); \quad\quad f_2=\frac { 1 }{ \pi \theta \left[ { 1+\left( \frac { x-\eta }{ \theta } \right) }^{ 2 } \right] },\quad \theta=1,\eta=0,x\in (1,\infty )\]

Problem 5.14

Obtain a Monte Carlo estimate of \(\int _{ 1 }^{ \infty }{ \frac { { x }^{ 2 } }{ \sqrt { 2\pi } } { e }^{ { -{ x }^{ 2 } }/{ 2 } } } dx\) by importance sampling.

g <- function(x) { (x^2 * exp(-x^2 / 2)) / sqrt(2 *pi) * (x > 1)}
integrate(g, 1, Inf)
## 0.400626 with absolute error < 5.7e-07
f1 <- function(x) { exp(-x) }
f2 <- function(x) { (pi * (1 + x^2))^(-1) * (x >= 1) }

m <- 10^7
x1 <- rexp(m)
x2 <- rcauchy(m)
x2[which(x2 < 1)] <- 1 # to catch overflow errors in g(x)

fg <- cbind(g(x1) / f1(x1), g(x2) / f2(x2))

theta.hat <- se <- numeric(2)
theta.hat <- c(mean(fg[,1]), mean(fg[,2]))
se <- c(sd(fg[,1]), sd(fg[,2]))
rbind(theta.hat, se)
##               [,1]      [,2]
## theta.hat 0.400505 0.4008194
## se        0.585777 0.9587950

Assignment Problems

Complete any one of the Kelton problems 6, 7, and 8.

Problem 7.9.6

A small free clinic has a single doctor who sees patients between 8:00 a.m. and noon. She spends 6 to 14 minutes (an average of 10) with each patient and so can theoretically see 6 patients per hour or a total of 24 each day. Staff is currently scheduling patient arrivals every 10 minutes, but have found that patients arrive as much as 15 minutes early or 30 minutes late, thus causing disruption to other patients. Worse, about 10 percent of patients never show up at all causing the doctor’s time to be wasted and wasting an appointment that could have been used by other sick patients. Staff would like to evaluate alternative scheduling strategies such as scheduling 2-3 arrivals every 20 minutes with the main objective of maximizing the doctors utilization (e.g., she would like to really see 24 patients each day if possible). Assume that the doctor will stay until all scheduled patients are seen, but she is very unhappy if she needs to stay past 12:30. Measure system performance mainly on the number of patients actually seen, but also consider the average patient waiting time, and how late the doctor must typically stay to see all scheduled patients.

Initial Setup

Using the AppointmentArrivals.spfx SimBit example:

  • Add Data Table under Data tab with:
    • Name of DailyAppointments.
    • Standard Property of Date-Time with:
      • Name of AppointmentTime.
      • Rows with data every 10 minutes from 8:00:00 AM to 11:50:00 PM.
  • Source with Entity Arrival Logic of:
    • Arrival Mode of Arrival Table
    • Arrival Time Property of DailyAppointments.AppointmentTime.
    • Arrival Time Deviation of Random.Uniform( -15 , 30 ) in Minutes.
    • Arrival No-Show Probability of 0.10.
  • Server with Processing Time of Random.Uniform( 6 , 14 ).
  • Experiment with Add Response of:
    • PatientsSeen using Doctor.Processing.NumberExited.
    • TimeInSystem using DefaultEntity.Population.TimeInSystem.Average.
    • DoctorUtilization using Doctor.Capacity.ScheduledUtilization.
    • HoursWorked using Doctor.ResourceState.TotalTime(1).
Model Data Table

Results

Alternatives

Starting with the above model, the following changes were made:

  • Add Data Table under Data tab with:
    • Name of AlternativeScheduling.
    • Standard Property of Date-Time with:
      • Name of EveryTwentyMinutes.
      • Rows with data every 20 minutes from 8:00:00 AM to 11:40:00 PM.
  • Source with Entity Arrival Logic of:
    • Arrival Mode of ArrivalTIme
    • Arrival Time Property of EntitiesPerArrival.
Results

Problem 7.9.7

Of great concern to emergency-department staff are those patients who on arrival are classified Severe but whose condition before receiving treatment deteriorates to Urgent, requiring immediate attention and stabilization in the Trauma Rooms. If such deterioration is noticed during the examination stage in 10% of the Severe patients, adjust Model 7-2 to reflect this and compare the impact on waiting times and throughput. What assumptions have you made about immediate attention?

Added new Entity of SevereUrgent. Added Path from ExamRooms output to TraumaRooms input.

  • Changed PatientMix for Severe in Data tab from \(24\) to \(24(1-0.1)=21.6\).
  • Changed Priority for Urgent in Data tab to \(5\).
  • Added PatientMix for SevereUrgent in Data tab of \(24(0.1)=2.4\).
  • Added Priority for SevereUrgent in Data tab of \(4\).
  • Added Treatments to SevereUrgent in Data tab of:
    • Sequence: Input@SignIn; Service Time: 1.
    • Sequence: Input@Registration; Service Time: 2.
    • Sequence: Input@ExamRooms; Service Time: Random.Triangular(15,20,25).
    • Sequence: Input@TraumaRooms; Service Time: Random.Triangular(15,25,35).
    • Sequence: Input@TreatmentRooms; Service Time: Random.Triangular(15,45,90).
    • Sequence: Input@TraumaExit; Service Time: 0.0.
Model Data Tab
Old Model New Model

Problem 7.9.8

Emergency-department staff work 8.5 hour shifts according to an “ideal” schedule that includes a one-hour meal break each shift - preferably around mid-shift. During meal breaks, services in each section are maintained. One of “Registration” staff members will cover for the “Sign In” officer during her meal break. The minimum staffing capacity of each area is given in [the below table]. When a new shift arrives there is a 30 minute briefing period" during which the two shifts team together to share information and handover. Develop a defensible staffing schedule and update example Model 7-2. Re-estimate staff utilization. What assumptions do you make concerning “handover time?”

Service Area Sign-In Registration Examination Treatment Rooms Trauma Rooms
Minimum Capacity 1 2 4 4 1

Three work schedules per 24-hour period were created with 1 hour overlap between each work schedule. The first 15-minute overlap is for the arriving worker(s) to don their uniform. The next 30-minute overlap was for information exchange to facilitate changeover. The last 15-minute overlap is for the departing worker(s) to change out of their uniform. The model was modified by changing Capacity Type from Fixed to Work Schedule with:

  • Min1Schedule based on DayPattern1 for SignIn Server.
  • Min2Schedule based on DayPattern2 for Registration Server.
  • Min4Schedule based on DayPattern3 for Examination Server.
  • Min4Schedule based on DayPattern3 for TreatmentRooms Server.
  • Min1Schedule based on DayPattern1 for TraumaRooms Server.
DayPattern1 DayPattern2 DayPattern3
Old Model New Model

References

http://s1.daumcdn.net/editor/fp/service_nc/pencil/Pencil_chromestore.html

https://en.wikipedia.org/wiki/Rayleigh_distribution#Generating_random_variates

Simio and Simulation: Modeling, Analysis, Applications 3d Ed. by W. David Kelton, Jeffrey S. Smith and David T. Sturrock with Simio software.

Discrete-Event Systems Simulation, 5th Edition (2010), by Jerry Banks, John S. Carlson, Barry L. Nelson,and David M. Nicol.