Simulasi Monte carlo di R

Saya mencoba menyelesaikan latihan berikut ini tetapi saya sebenarnya tidak tahu bagaimana memulai melakukan ini. Saya telah menemukan beberapa kode dalam buku saya yang terlihat seperti itu tetapi ini adalah latihan yang sama sekali berbeda dan saya tidak tahu bagaimana menghubungkannya satu sama lain. Bagaimana saya bisa mulai mensimulasikan kedatangan dan bagaimana saya tahu kapan mereka selesai? Saya tahu cara menyimpannya dan menghitung a, b, c, d sesuai dengan itu. Tapi saya tidak tahu bagaimana sebenarnya saya perlu mensimulasikan simulasi monte carlo. Bisakah seseorang membantu saya memulai? Saya tahu ini bukan tempat di mana pertanyaan Anda dijawab untuk Anda tetapi hanya diselesaikan sebagai gantinya. Tetapi masalahnya adalah saya tidak tahu bagaimana memulainya.

Meja bantuan dukungan TI mewakili sistem antrian dengan lima asisten menerima panggilan dari pelanggan. Panggilan terjadi sesuai dengan proses Poisson dengan tingkat rata-rata satu panggilan setiap 45 detik. Waktu layanan untuk asisten 1, 2, 3, 4, dan 5 adalah semua variabel acak Eksponensial dengan parameter λ1 = 0,1, λ2 = 0,2, λ3 = 0,3, λ4 = 0,4, dan λ5 = 0,5 mnt − 1, masing-masing ( jth help desk assistant memiliki λk = k / 10 min − 1). Selain pelanggan yang sedang dibantu, hingga sepuluh pelanggan lain dapat ditahan. Pada saat kapasitas ini tercapai, penelepon baru menerima sinyal sibuk. Gunakan metode Monte Carlo untuk memperkirakan karakteristik kinerja berikut,

(a) sebagian kecil dari pelanggan yang menerima sinyal sibuk;

(B) waktu respon yang diharapkan;

(c) waktu tunggu rata-rata;

(d) porsi pelanggan yang dilayani oleh setiap asisten help desk;

EDIT: yang saya miliki sejauh ini (tidak banyak):

pa = 1/45sec-1

jobs = rep(1,5); onHold = rep(1,10);

jobsIndex = 0;

onHoldIndex = 0;

u = runif(1)
for (i in 1:1000) {

    if(u  <= pa){ # new arrival

        if(jobsIndex < 5) # assistant is free, #give job to assistant

            jobsIndex++;

        else #add to onHold array

            onHoldIndex++;
    }
}

Ini adalah salah satu jenis simulasi yang paling instruktif dan menyenangkan untuk dilakukan: Anda membuat agen independen di komputer, membiarkan mereka berinteraksi, melacak apa yang mereka lakukan, dan mempelajari apa yang terjadi. Ini adalah cara yang luar biasa untuk belajar tentang sistem yang kompleks, terutama (tetapi tidak terbatas pada) yang tidak dapat dipahami dengan analisis matematika murni.

Cara terbaik untuk membangun simulasi semacam itu adalah dengan desain top-down.

Pada tingkat yang paling tinggi , kode tersebut akan terlihat seperti

initialize(...)
while (process(get.next.event())) {}

(Ini dan semua contoh berikutnya adalah kode yang dapat dieksekusi R , bukan hanya kode semu.) Loop adalah simulasi yang didorong oleh peristiwa : get.next.event()menemukan "peristiwa" yang menarik dan meneruskan deskripsi tentang hal itu process, yang melakukan sesuatu dengannya (termasuk mencatat semua informasi tentang itu). Ia kembali TRUEselama hal-hal berjalan dengan baik; setelah mengidentifikasi kesalahan atau akhir simulasi, ia kembali FALSE, mengakhiri loop.

Jika kita membayangkan implementasi fisik antrian ini, seperti orang yang menunggu surat nikah di New York City atau untuk surat izin mengemudi atau tiket kereta api hampir di mana saja, kita memikirkan dua jenis agen: pelanggan dan "asisten" (atau server) . Pelanggan mengumumkan diri mereka dengan muncul; asisten mengumumkan ketersediaan mereka dengan menyalakan lampu atau menandatangani atau membuka jendela. Ini adalah dua jenis acara untuk diproses.

Lingkungan ideal untuk simulasi semacam itu adalah lingkungan berorientasi objek yang benar di mana objek bisa berubah : mereka dapat mengubah keadaan untuk merespons secara independen terhadap hal-hal di sekitar mereka. Rbenar-benar mengerikan untuk ini (bahkan Fortran akan lebih baik!). Namun, kita masih bisa menggunakannya jika kita berhati-hati. Kuncinya adalah menjaga semua informasi dalam satu set umum struktur data yang dapat diakses (dan dimodifikasi) oleh banyak prosedur yang terpisah dan saling berinteraksi. Saya akan mengadopsi konvensi menggunakan nama variabel DALAM SEMUA CAPS untuk data tersebut.

Level selanjutnya dari desain top-down adalah kode process. Itu menanggapi deskriptor acara tunggal e:

process <- function(e) {
  if (is.null(e)) return(FALSE)
  if (e$type == "Customer") {
    i <- find.assistant(e$time)
    if (is.null(i)) put.on.hold(e$x, e$time) else serve(i, e$x, e$time)
  } else {
    release.hold(e$time)
  }
  return(TRUE)
}

Itu harus menanggapi acara nol ketika get.next.eventtidak ada acara untuk dilaporkan. Jika tidak, processterapkan "aturan bisnis" sistem. Praktis menulis sendiri dari uraian dalam pertanyaan. Cara kerjanya harus memerlukan sedikit komentar, kecuali untuk menunjukkan bahwa pada akhirnya kita perlu kode subrutin put.on.holddan release.hold(menerapkan antrian pelanggan-memegang) dan serve(menerapkan interaksi pelanggan-asisten).

Apa itu "acara"? Itu harus berisi informasi tentang siapa yang bertindak, tindakan apa yang mereka ambil, dan kapan itu terjadi. Karena itu kode saya menggunakan daftar yang berisi ketiga jenis informasi ini. Namun, get.next.eventhanya perlu memeriksa waktu. Ini bertanggung jawab hanya untuk menjaga antrian acara di mana

1. Setiap peristiwa dapat dimasukkan ke dalam antrian ketika diterima dan

2. Peristiwa paling awal dalam antrian dapat dengan mudah diekstraksi dan diteruskan ke pemanggil.

Implementasi terbaik dari antrian prioritas ini akan menjadi tumpukan, tapi itu terlalu cerewet R. Mengikuti saran dalam The Art of R Programming karya Norman Matloff (yang menawarkan simulator antrian yang lebih fleksibel, abstrak, tetapi terbatas), saya telah menggunakan kerangka data untuk mengadakan acara dan hanya mencari waktu minimum di antara catatan-catatannya.

get.next.event <- function() {
  if (length(EVENTS$time) <= 0) new.customer()               # Wait for a customer$
  if (length(EVENTS$time) <= 0) return(NULL)                 # Nothing's going on!$
  if (min(EVENTS$time) > next.customer.time()) new.customer()# See text
  i <- which.min(EVENTS$time)
  e <- EVENTS[i, ]; EVENTS <<- EVENTS[-i, ]
  return (e)
}

Ada banyak cara ini bisa dikodekan. Versi terakhir yang diperlihatkan di sini mencerminkan pilihan yang saya buat dalam mengodekan bagaimana processbereaksi terhadap peristiwa "Asisten" dan cara new.customerkerjanya: get.next.eventhanya mengeluarkan pelanggan dari antrean penangguhan, lalu duduk kembali dan menunggu acara lain. Terkadang perlu mencari pelanggan baru dengan dua cara: pertama, untuk melihat apakah ada yang menunggu di pintu (seolah-olah) dan kedua, apakah ada yang masuk ketika kita tidak melihat.

Jelas, new.customerdan next.customer.timeini rutinitas penting , jadi mari kita urus mereka selanjutnya.

new.customer <- function() {  
  if (CUSTOMER.COUNT < dim(CUSTOMERS)[2]) {
    CUSTOMER.COUNT <<- CUSTOMER.COUNT + 1
    insert.event(CUSTOMER.COUNT, "Customer", 
                 CUSTOMERS["Arrived", CUSTOMER.COUNT])
  }
  return(CUSTOMER.COUNT)
}
next.customer.time <- function() {
  if (CUSTOMER.COUNT < dim(CUSTOMERS)[2]) {
    x <- CUSTOMERS["Arrived", CUSTOMER.COUNT]
  } else {x <- Inf}
  return(x) # Time when the next customer will arrive
}

CUSTOMERSadalah array 2D, dengan data untuk setiap pelanggan di kolom. Ini memiliki empat baris (bertindak sebagai bidang) yang menggambarkan pelanggan dan mencatat pengalaman mereka selama simulasi : "Tiba", "Dilayani", "Durasi", dan "Asisten" (pengidentifikasi numerik positif dari asisten, jika ada, yang melayani mereka, dan sebaliknya -1untuk sinyal sibuk). Dalam simulasi yang sangat fleksibel, kolom-kolom ini akan dihasilkan secara dinamis, tetapi karena Rsuka bekerja, mudah untuk menghasilkan semua pelanggan pada awalnya, dalam satu matriks besar, dengan waktu kedatangan mereka sudah dihasilkan secara acak. next.customer.timedapat mengintip kolom berikutnya dari matriks ini untuk melihat siapa yang akan datang berikutnya. Variabel globalCUSTOMER.COUNTmenunjukkan pelanggan terakhir yang tiba. Pelanggan dikelola dengan sangat sederhana melalui penunjuk ini, memajukannya untuk mendapatkan pelanggan baru dan melihat melampaui itu (tanpa memajukan) untuk mengintip pelanggan berikutnya.

serve mengimplementasikan aturan bisnis dalam simulasi.

serve <- function(i, x, time.now) {
  #
  # Serve customer `x` with assistant `i`.
  #
  a <- ASSISTANTS[i, ]
  r <- rexp(1, a$rate)                       # Simulate the duration of service
  r <- round(r, 2)                           # (Make simple numbers)
  ASSISTANTS[i, ]$available <<- time.now + r # Update availability
  #
  # Log this successful service event for later analysis.
  #
  CUSTOMERS["Assistant", x] <<- i
  CUSTOMERS["Served", x] <<- time.now
  CUSTOMERS["Duration", x] <<- r
  #
  # Queue the moment the assistant becomes free, so they can check for
  # any customers on hold.
  #
  insert.event(i, "Assistant", time.now + r)
  if (VERBOSE) cat(time.now, ": Assistant", i, "is now serving customer", 
                   x, "until", time.now + r, "\n")
  return (TRUE)
}

Ini mudah. ASSISTANTSadalah kerangka data dengan dua bidang: capabilities(memberikan tingkat layanan mereka) dan available, yang menandai waktu berikutnya di mana asisten akan bebas. Pelanggan dilayani dengan menghasilkan durasi layanan acak sesuai dengan kemampuan asisten, memperbarui waktu ketika asisten berikutnya tersedia, dan mencatat interval layanan dalam CUSTOMERSstruktur data. The VERBOSEbendera berguna untuk pengujian dan debugging: ketika benar, itu memancarkan aliran kalimat bahasa Inggris yang menggambarkan poin pengolahan kunci.

Bagaimana asisten ditugaskan kepada pelanggan adalah penting dan menarik. Seseorang dapat membayangkan beberapa prosedur: penugasan secara acak, dengan beberapa pemesanan tetap, atau menurut siapa yang telah bebas waktu paling lama (atau terpendek). Banyak dari ini diilustrasikan dalam kode komentar:

find.assistant <- function(time.now) {
  j <- which(ASSISTANTS$available <= time.now)
  #if (length(j) > 0) {
  #  i <- j[ceiling(runif(1) * length(j))]
  #} else i <- NULL                                    # Random selection
  #if (length(j) > 0) i <- j[1] else i <- NULL         # Pick first assistant
  #if (length(j) > 0) i <- j[length(j)] else i <- NULL # Pick last assistant
  if (length(j) > 0) {
    i <- j[which.min(ASSISTANTS[j, ]$available)]
  } else i <- NULL                                     # Pick most-rested assistant
  return (i)
}

Sisa simulasi sebenarnya hanya latihan rutin dalam membujuk Runtuk menerapkan struktur data standar, terutama buffer melingkar untuk antrian yang ditahan. Karena Anda tidak ingin mengamuk dengan global, saya menempatkan semua ini dalam satu prosedur sim. Argumennya menggambarkan masalah: jumlah pelanggan yang disimulasikan ( n.events), tingkat kedatangan pelanggan, kemampuan asisten, dan ukuran antrian tahan (yang dapat diatur ke nol untuk menghilangkan antrian sama sekali).

r <- sim(n.events=250, arrival.rate=60/45, capabilities=1:5/10, hold.queue.size=10)

Ini mengembalikan daftar struktur data yang dipertahankan selama simulasi; salah satu yang paling menarik adalah CUSTOMERSarray. Rmembuatnya cukup mudah untuk merencanakan informasi penting dalam array ini dengan cara yang menarik. Berikut adalah satu output yang menunjukkan pelanggan terakhir dalam simulasi yang lebih lama dari pelanggan.$50$$250$

Pengalaman setiap pelanggan diplot sebagai garis waktu horizontal, dengan simbol lingkaran pada saat kedatangan, garis hitam solid untuk setiap penantian yang menunggu, dan garis berwarna untuk durasi interaksi mereka dengan asisten (warna dan jenis garis) bedakan di antara asisten) Di bawah plot Pelanggan ini adalah salah satu yang menunjukkan pengalaman asisten, menandai waktu ketika mereka dan tidak terlibat dengan pelanggan. Titik akhir setiap interval aktivitas dibatasi oleh batang vertikal.

Saat dijalankan verbose=TRUE, output teks simulasi terlihat seperti ini:

...
160.71 : Customer 211 put on hold at position 1 
161.88 : Customer 212 put on hold at position 2 
161.91 : Assistant 3 is now serving customer 213 until 163.24 
161.91 : Customer 211 put on hold at position 2 
162.68 : Assistant 4 is now serving customer 212 until 164.79 
162.71 : Assistant 5 is now serving customer 211 until 162.9 
163.51 : Assistant 5 is now serving customer 214 until 164.05 
...

(Angka di sebelah kiri adalah waktu setiap pesan dipancarkan.) Anda dapat mencocokkan deskripsi ini dengan bagian plot Pelanggan yang terbentang di antara waktu dan .$160$$165$

Kami dapat mempelajari pengalaman pelanggan yang ditahan dengan memplot durasi yang ditahan oleh pengidentifikasi pelanggan, menggunakan simbol khusus (merah) untuk menunjukkan kepada pelanggan yang menerima sinyal sibuk.

(Bukankah semua plot ini menghasilkan dasbor real-time yang hebat bagi siapa pun yang mengelola antrian layanan ini!)

Sangat menarik untuk membandingkan plot dan statistik yang Anda dapatkan dengan memvariasikan parameter yang dilewati sim. Apa yang terjadi ketika pelanggan datang terlalu cepat untuk diproses? Apa yang terjadi ketika antrian terus dibuat lebih kecil atau dihilangkan? Apa yang berubah ketika asisten dipilih dengan cara berbeda? Bagaimana angka dan kemampuan asisten mempengaruhi pengalaman pelanggan? Apa poin penting di mana beberapa pelanggan mulai ditolak atau mulai ditunda untuk waktu yang lama?

Biasanya, untuk pertanyaan belajar mandiri yang jelas seperti ini, kami akan berhenti di sini dan meninggalkan detail yang tersisa sebagai latihan. Namun, saya tidak ingin mengecewakan pembaca yang mungkin sudah sejauh ini dan tertarik untuk mencoba ini sendiri (dan mungkin memodifikasi dan membangunnya untuk tujuan lain), jadi di bawah ini adalah kode kerja lengkap.

( Pemrosesan di situs ini akan mengacaukan lekukan pada setiap baris yang mengandung simbol , tetapi lekukan yang dapat dibaca harus dikembalikan ketika kode disisipkan ke dalam file teks.)$\TeX$$\$$

sim <- function(n.events, verbose=FALSE, ...) {
  #
  # Simulate service for `n.events` customers.
  #
  # Variables global to this simulation (but local to the function):
  #
  VERBOSE <- verbose         # When TRUE, issues informative message
  ASSISTANTS <- list()       # List of assistant data structures
  CUSTOMERS <- numeric(0)    # Array of customers that arrived
  CUSTOMER.COUNT <- 0        # Number of customers processed
  EVENTS <- list()           # Dynamic event queue   
  HOLD <- list()             # Customer on-hold queue
  #............................................................................#
  #
  # Start.
  #
  initialize <- function(arrival.rate, capabilities, hold.queue.size) {
    #
    # Create common data structures.
    #
    ASSISTANTS <<- data.frame(rate=capabilities,     # Service rate
                              available=0            # Next available time
    )
    CUSTOMERS <<- matrix(NA, nrow=4, ncol=n.events, 
                         dimnames=list(c("Arrived",  # Time arrived
                                         "Served",   # Time served
                                         "Duration", # Duration of service
                                         "Assistant" # Assistant id
                         )))
    EVENTS <<- data.frame(x=integer(0),              # Assistant or customer id
                          type=character(0),         # Assistant or customer
                          time=numeric(0)            # Start of event
    )
    HOLD <<- list(first=1,                           # Index of first in queue
                  last=1,                            # Next available slot
                  customers=rep(NA, hold.queue.size+1))
    #
    # Generate all customer arrival times in advance.
    #
    CUSTOMERS["Arrived", ] <<- cumsum(round(rexp(n.events, arrival.rate), 2))
    CUSTOMER.COUNT <<- 0
    if (VERBOSE) cat("Started.\n")
    return(TRUE)
  }
  #............................................................................#
  #
  # Dispatching.
  #
  # Argument `e` represents an event, consisting of an assistant/customer 
  # identifier `x`, an event type `type`, and its time of occurrence `time`.
  #
  # Depending on the event, a customer is either served or an attempt is made
  # to put them on hold.
  #
  # Returns TRUE until no more events occur.
  #
  process <- function(e) {
    if (is.null(e)) return(FALSE)
    if (e$type == "Customer") {
      i <- find.assistant(e$time)
      if (is.null(i)) put.on.hold(e$x, e$time) else serve(i, e$x, e$time)
    } else {
      release.hold(e$time)
    }
    return(TRUE)
  }#$
  #............................................................................#
  #
  # Event queuing.
  #
  get.next.event <- function() {
    if (length(EVENTS$time) <= 0) new.customer()
    if (length(EVENTS$time) <= 0) return(NULL)
    if (min(EVENTS$time) > next.customer.time()) new.customer()
    i <- which.min(EVENTS$time)
    e <- EVENTS[i, ]; EVENTS <<- EVENTS[-i, ]
    return (e)
  }
  insert.event <- function(x, type, time.occurs) {
    EVENTS <<- rbind(EVENTS, data.frame(x=x, type=type, time=time.occurs))
    return (NULL)
  }
  # 
  # Customer arrivals (called by `get.next.event`).
  #
  # Updates the customers pointer `CUSTOMER.COUNT` and returns the customer
  # it newly points to.
  #
  new.customer <- function() {  
    if (CUSTOMER.COUNT < dim(CUSTOMERS)[2]) {
      CUSTOMER.COUNT <<- CUSTOMER.COUNT + 1
      insert.event(CUSTOMER.COUNT, "Customer", 
                   CUSTOMERS["Arrived", CUSTOMER.COUNT])
    }
    return(CUSTOMER.COUNT)
  }
  next.customer.time <- function() {
    if (CUSTOMER.COUNT < dim(CUSTOMERS)[2]) {
      x <- CUSTOMERS["Arrived", CUSTOMER.COUNT]
    } else {x <- Inf}
    return(x) # Time when the next customer will arrive
  }
  #............................................................................#
  #
  # Service.
  #
  find.assistant <- function(time.now) {
    #
    # Select among available assistants.
    #
    j <- which(ASSISTANTS$available <= time.now) 
    #if (length(j) > 0) {
    #  i <- j[ceiling(runif(1) * length(j))]
    #} else i <- NULL                                    # Random selection
    #if (length(j) > 0) i <- j[1] else i <- NULL         # Pick first assistant
    #if (length(j) > 0) i <- j[length(j)] else i <- NULL # Pick last assistant
    if (length(j) > 0) {
      i <- j[which.min(ASSISTANTS[j, ]$available)]
    } else i <- NULL # Pick most-rested assistant
    return (i)
  }#$
  serve <- function(i, x, time.now) {
    #
    # Serve customer `x` with assistant `i`.
    #
    a <- ASSISTANTS[i, ]
    r <- rexp(1, a$rate)                       # Simulate the duration of service
    r <- round(r, 2)                           # (Make simple numbers)
    ASSISTANTS[i, ]$available <<- time.now + r # Update availability
    #
    # Log this successful service event for later analysis.
    #
    CUSTOMERS["Assistant", x] <<- i
    CUSTOMERS["Served", x] <<- time.now
    CUSTOMERS["Duration", x] <<- r
    #
    # Queue the moment the assistant becomes free, so they can check for
    # any customers on hold.
    #
    insert.event(i, "Assistant", time.now + r)
    if (VERBOSE) cat(time.now, ": Assistant", i, "is now serving customer", 
                     x, "until", time.now + r, "\n")
    return (TRUE)
  }
  #............................................................................#
  #
  # The on-hold queue.
  #
  # This is a cicular buffer implemented by an array and two pointers,
  # one to its head and the other to the next available slot.
  #
  put.on.hold <- function(x, time.now) {
    #
    # Try to put customer `x` on hold.
    #
    if (length(HOLD$customers) < 1 || 
          (HOLD$first - HOLD$last %% length(HOLD$customers) == 1)) {
      # Hold queue is full, alas.  Log this occurrence for later analysis.
      CUSTOMERS["Assistant", x] <<- -1 # Busy signal
      if (VERBOSE) cat(time.now, ": Customer", x, "got a busy signal.\n")
      return(FALSE)
    }
    #
    # Add the customer to the hold queue.
    #
    HOLD$customers[HOLD$last] <<- x
    HOLD$last <<- HOLD$last %% length(HOLD$customers) + 1
    if (VERBOSE) cat(time.now, ": Customer", x, "put on hold at position", 
                 (HOLD$last - HOLD$first - 1) %% length(HOLD$customers) + 1, "\n")
    return (TRUE)
  }
  release.hold <- function(time.now) {
    #
    # Pick up the next customer from the hold queue and place them into
    # the event queue.
    #
    if (HOLD$first != HOLD$last) {
      x <- HOLD$customers[HOLD$first]   # Take the first customer
      HOLD$customers[HOLD$first] <<- NA # Update the hold queue
      HOLD$first <<- HOLD$first %% length(HOLD$customers) + 1
      insert.event(x, "Customer", time.now)
    }
  }$
  #............................................................................#
  #
  # Summaries.
  #
  # The CUSTOMERS array contains full information about the customer experiences:
  # when they arrived, when they were served, how long the service took, and
  # which assistant served them.
  #
  summarize <- function() return (list(c=CUSTOMERS, a=ASSISTANTS, e=EVENTS,
                                       h=HOLD))
  #............................................................................#
  #
  # The main event loop.
  #
  initialize(...)
  while (process(get.next.event())) {}
  #
  # Return the results.
  #
  return (summarize())
}
#------------------------------------------------------------------------------#
#
# Specify and run a simulation.
#
set.seed(17)
n.skip <- 200  # Number of initial events to skip in subsequent summaries
system.time({
  r <- sim(n.events=50+n.skip, verbose=TRUE, 
           arrival.rate=60/45, capabilities=1:5/10, hold.queue.size=10)
})
#------------------------------------------------------------------------------#
#
# Post processing.
#
# Skip the initial phase before equilibrium.
#
results <- r$c
ids <- (n.skip+1):(dim(results)[2])
arrived <- results["Arrived", ]
served <- results["Served", ]
duration <- results["Duration", ]
assistant <- results["Assistant", ]
assistant[is.na(assistant)] <- 0   # Was on hold forever
ended <- served + duration
#
# A detailed plot of customer experiences.
#
n.events <- length(ids)
n.assistants <- max(assistant, na.rm=TRUE) 
colors <- rainbow(n.assistants + 2)
assistant.color <- colors[assistant + 2]
x.max <- max(results["Served", ids] + results["Duration", ids], na.rm=TRUE)
x.min <- max(min(results["Arrived", ids], na.rm=TRUE) - 2, 0)
#
# Lay out the graphics.
#
layout(matrix(c(1,1,2,2), 2, 2, byrow=TRUE), heights=c(2,1))
#
# Set up the customers plot.
#
plot(c(x.min, x.max), range(ids), type="n",
     xlab="Time", ylab="Customer Id", main="Customers")
#
# Place points at customer arrival times.
#
points(arrived[ids], ids, pch=21, bg=assistant.color[ids], col="#00000070")
#
# Show wait times on hold.
#
invisible(sapply(ids, function(i) {
  if (!is.na(served[i])) lines(x=c(arrived[i], served[i]), y=c(i,i))
}))
#
# More clearly show customers getting a busy signal.
#
ids.not.served <- ids[is.na(served[ids])]
ids.served <- ids[!is.na(served[ids])]
points(arrived[ids.not.served], ids.not.served, pch=4, cex=1.2)
#
# Show times of service, colored by assistant id.
#
invisible(sapply(ids.served, function(i) {
  lines(x=c(served[i], ended[i]), y=c(i,i), col=assistant.color[i], lty=assistant[i])
}))
#
# Plot the histories of the assistants.
#
plot(c(x.min, x.max), c(1, n.assistants)+c(-1,1)/2, type="n", bty="n",
     xlab="", ylab="Assistant Id", main="Assistants")
abline(h=1:n.assistants, col="#808080", lwd=1)
invisible(sapply(1:(dim(results)[2]), function(i) {
  a <- assistant[i]
  if (a > 0) {
    lines(x=c(served[i], ended[i]), y=c(a, a), lwd=3, col=colors[a+2])
    points(x=c(served[i], ended[i]), y=c(a, a), pch="|", col=colors[a+2])
  }
}))
#
# Plot the customer waiting statistics.
#
par(mfrow=c(1,1))
i <- is.na(served)
plot(served - arrived, xlab="Customer Id", ylab="Minutes",
     main="Service Wait Durations")
lines(served - arrived, col="Gray")
points(which(i), rep(0, sum(i)), pch=16, col="Red")
#
# Summary statistics.
#
mean(!is.na(served)) # Proportion of customers served
table(assistant)