Statistical Analysis of a Telephone Call Center: A Queueing

15 Pages·2005·0.43 MB·English

Statistical Analysis of a Telephone Call Center:
A Queueing-Science Perspective
Lawrence BROWN, Noah GANS, Avishai MANDELBAUM, Anat SAKOV,
Haipeng SHEN, Sergey ZELTYN, andLinda ZHAO
Acallcenterisaservicenetworkinwhichagentsprovidetelephone-based services.Customerswhoseektheseservicesaredelayedin
tele-queues. Thisarticlesummarizesananalysisofauniquerecordofcallcenteroperations.Thedatacompriseacompleteoperational
historyofasmallbankingcallcenter,callbycall,overafullyear.Takingtheperspectiveofqueueingtheory,wedecomposetheservice
processintothreefundamentalcomponents:arrivals,customerpatience, andservicedurations.Eachcomponentinvolvesdifferentbasic
mathematical structures andrequires adifferent styleofstatistical analysis.Someofthekeyempirical results aresketched, alongwith
descriptions ofthevariedtechniques required.Severalstatisticaltechniques aredevelopedforanalysisofthebasiccomponents.Oneof
thesetechniquesisatestthatapointprocessisaPoissonprocess.Anotherinvolvesestimationofthemeanfunctioninanonparametric
regressionwithlognormalerrors.Anewgraphicaltechniqueisintroducedfornonparametrichazardrateestimationwithcensoreddata.
ModelsaredevelopedandimplementedforforecastingofPoissonarrivalrates.Finally,thearticlesurveyshowthecharacteristicsdeduced
fromthestatisticalanalysesformthebuildingblocksfortheoreticallyinterestingandpracticallyusefulmathematicalmodelsforcallcenter
operations.
KEY WORDS: Abandonment; Arrivals; Call center; Censored data; Erlang-A; Erlang-C; Human patience; Inhomogeneous Poisson
process;Khintchine–Pollaczekformula;Lognormaldistribution;Multiserverqueue;PredictionofPoissonrates;Queue-
ingscience;Queueingtheory;Servicetime.
1. INTRODUCTION levelofdetail,havebeenscarce.Thedatathataretypicallycol-
lectedandusedinthecallcenterindustryaresimpleaverages
Telephone call centers are technology-intensive operations.
calculatedforthecallsthatarrivewithinfixedintervalsoftime,
Nevertheless, often 70% or more of their operating costs are
often 30 minutes. There is a lack of documented,comprehen-
devoted to human resources. Well-run call centers adhere to
sive, empirical research on call center performance that uses
asharplydefinedbalancebetweenagentefficiencyandservice
moredetaileddata.
quality; to do so, they use queueing-theoretic models. Inputs
Theimmediategoalofourstudyistofillthisgap.Inthisarti-
to these mathematicalmodelsare statistics concerningsystem
cle,wesummarizeacomprehensiveanalysisofoperationaldata
primitives, such as the number of agents working, the rate at
from a bank call center. The data span all 12 months of 1999
whichcallsarrive,thetimerequiredforacustomertobeserved,
andarecollectedatthelevelofindividualcalls.Ourdatasource
and the length of time customers are willing to wait on hold
consistsofmorethan1,200,000callsthatarrivedatthecenter
before they hang up the phone and abandon the queue. Out-
over the year. Of these, about 750,000 calls terminated in an
putsareperformancemeasures,suchasthedistributionoftime
interactive voice response unit (IVR or VRU), a type of an-
thatcustomerswait“onhold”andthefractionofcustomersthat
swering machine that allows customers to serve themselves.
abandonthequeuebeforebeingserved.Inpractice,thenumber
Theremaining450,000callersaskedtobeservedbyanagent;
of agentsworkingbecomesa controlparameter,which canbe
wehavearecordoftheevent-historyofeachofthesecalls.
increased or decreased to attain the desired efficiency–quality
This article is an important part of a larger effort to use
trade-off.
both theoretical and empirical tools to better characterize call
Estimatesofthese primitivesare neededto calibratequeue-
center operations and performance. It is an abridged version
ing models, and in many cases the models make distribu-
of the work of Brown et al. (2002a), which provided a more
tional assumptions concerning the primitives. In theory, the
complete treatment of the results reported here. Mandelbaum,
datarequiredtovalidateandproperlytunethesemodelsshould
Sakov, and Zeltyn (2000)presented a comprehensivedescrip-
be readily available, because computers track and control the
tionofourcall-by-calldatabase.Gans,Koole,andMandelbaum
minutestdetailsofeverycall’sprogressthroughthesystem.Itis
(2003) reviewed queueing and related models of call centers,
thussurprisingthatoperationaldata,collectedatanappropriate
andMandelbaum(2001)providedanextensivebibliography.
1.1 QueueingModelsofCallCenters
Lawrence Brown is Professor, Department of Statistics, The Wharton
School,UniversityofPennsylvania,Philadelphia,PA19104(E-mail:lbrown@
The simplest and most widely used queueing model in call
wharton.upenn.edu).NoahGansisAssociateProfessor,DepartmentofOpera-
tionsandInformationManagement,TheWhartonSchool,UniversityofPenn- centersistheso-calledM/M/Nsystem,sometimesreferredto
sylvania,Philadelphia,PA19104(E-mail:[email protected]).Avishai asErlang-C(Erlang1911,1917).TheM/M/Nmodelisquite
MandelbaumisProfessor,FacultyofIndustrialEngineeringandManagement,
restrictive. It assumes, among other things, a steady-state en-
Technion,Haifa,Israel(E-mail:[email protected]). AnatSakovisPost-
doctoralFellow,Tel-AvivUniversity,Tel-Aviv,Israel(E-mail:[email protected]. vironmentinwhicharrivalsconformtoaPoissonprocess,ser-
ac.il).HaipengShenisAssistantProfessor,DepartmentofStatistics, Univer- vicedurationsareexponentiallydistributed,andcustomersand
sityofNorthCarolina,Durham,NC27599(E-mail:[email protected]).
serversarestatisticallyidenticalandactindependentlyofeach
SergeyZeltynisPh.D.Candidate,FacultyofIndustrialEngineeringandMan-
agement, Technion, Haifa, Israel (E-mail: [email protected]). Linda other. It does not acknowledge,among other things, customer
Zhao is Associate Professor, Department of Statistics, The Wharton School,
UniversityofPennsylvania,Philadelphia,PA19104(E-mail:lzhao@wharton.
upenn.edu). This work was supported by National Science Foundation ©2005AmericanStatisticalAssociation
DMS-99-71751 and DMS-99-71848, the Sloane Foundation, Israeli Science JournaloftheAmericanStatisticalAssociation
Foundationgrants388/99and126/02,theWhartonFinancialInstitutionsCen- March2005,Vol.100,No.469,ApplicationsandCaseStudies
ter,andTechnionfundsforthepromotionofresearchandsponsoredresearch. DOI10.1198/016214504000001808
36
Brownetal.:StatisticalAnalysisofaTelephoneCallCenter 37
impatience and abandonment behavior, time-dependent para- whoiscapableofperformingthedesiredservice,thenthecus-
meters, customers’ heterogeneity, or servers’ skill levels. An tomer and the agent are matched to start service immediately.
essentialtaskofcontemporaryqueueingtheoristsistodevelop Otherwise,thecustomerjoinsthetele-queue.
modelsthataccountfortheseeffects. Customersinthetele-queuearenominallyservedonafirst-
Queueing science seeks to determine which of these ef- come, first-served (FCFS) basis, and customers’ positions in
fects is most important for modeling real-life situations. For queuearedistinguishedbythetimeswhentheyarrive.Inprac-
example, Garnett, Mandelbaum, and Reiman (2002) devel- tice, the call center operates a system with two priorities—
opedbothexactandapproximateexpressionsforM/M/N+M highandlow—andmoveshigh-prioritycustomersupinqueue
(also called Erlang-A) systems, which explicitly model cus- by subtracting 1.5 minutes from their actual arrival times.
tomer patience (time to abandonment) as being exponentially Mandelbaum et al. (2000) compared the behavior of the two
distributed. Empirical analysis can help us judge how well prioritygroupsofcustomers.
the Erlang-C and Erlang-A models predict customer delays, Whilewaiting,eachcustomerperiodicallyreceivesinforma-
whetherornottheirunderlyingassumptionsaremet. tion on his or her progressin the queue.More specifically,he
orsheistoldtheamountoftimethatthefirstpersonin queue
1.2 StructureoftheArticle
hasbeenwaiting,aswellashisorherapproximatelocationin
The article is structured as follows. Section 2 describes the thequeue.Theannouncementisreplayedevery60secondsor
callcenterunderstudyand its database.Eachof Sections3–5 so,withmusic,news,orcommercialsintertwined.
is dedicated to the statistical analysis of one of the stochastic Ineachofthe12monthsof1999,roughly100,000–120,000
primitivesofthequeueingsystem:Section3addressescallar- calls arrived to the system, with 65,000–85,000 of these ter-
rivals;Section4,servicedurations;andSection5,tele-queueing minating in the VRU. The remaining30,000–40,000calls per
andcustomerpatience.Section5alsoanalyzescustomerwait- monthinvolvedcallerswhoexitedtheVRUindicatingadesire
ingtimes, a performancemeasuredeeplyintertwinedwith the to speak to an agent. These calls are the focus of our study.
abandonmentprimitive. About80%ofthoserequestingservicewereinfactserved,and
A synthesis of the primitive building blocks is typically about20%wereabandonedbeforebeingserved.
needed for operational understanding. Toward this end, Sec- Each call that proceeds past the VRU can be thought of as
tion6discussespredictionofthearriving“workload,”whichis passingthroughuptothreestages,eachofwhichgeneratesdis-
essentialinpracticeforsettingsuitableservicestaffinglevels. tinctdata. Thefirst of these is the arrivalstage, which is trig-
Once each of the primitives has been analyzed, one can gered by the call’s exit from the VRU and generates a record
also attempt to use existing queueing theory,or modifications ofanarrivaltime.Ifnoappropriateserverisavailable,thenthe
thereof,to describe certain featuresof the holistic behaviorof callentersthequeueingstage.Threepiecesofdataarerecorded
the system. Section 7 concludes with analyses of this type. foreachcallthatqueues:thetimeitenteredthequeue,thetime
We validate some classical theoretical results from queueing itexitedthequeue,andthemannerinwhichitexitedthequeue,
theoryandrefuteothers. by being served or abandoning. In the last stage, service, the
Finally, we note that many statistical tests are considered data recordedare the starting and endingtimesof the service.
throughoutthearticle,whichraisestheproblemofmultiplicity Note that calls that are served immediately skip the queueing
(Benjaminiand Hochberg1995).When data fromcall centers stage,andcallsthatareabandonedneverentertheservicestage.
areanalyzedinsupportofoperationaldecisions,themultiplic- Inadditiontothesetimestamps,eachcallrecordinourdata-
ityproblemmustbeaddressed. base includes a categorical description of the type of service
requested. The main call types are regular (PS in the data-
2. THECALLCENTEROFBANKANONYMOUS
base), stock transaction (NE), new/potential customer (NW),
The source of our data (Call Center Data 2002) is a small and Internet assistance (IN). Mandelbaum et al. (2000) de-
call centerfor oneof Israel’s banks.This centerprovidessev- scribedtheprocessofcollectingandcleaningthedataandpro-
eral types of basic services, as well as others, including stock videdadditionaldescriptiveanalysisofthedata.
trading and technical support, for users of the bank’s Internet Over the year, two importantoperationalchangesoccurred.
site. On weekdays (Sunday–Thursday in Israel) the center is First,inJanuary–July,allcallswereservedbythesamegroup
open from 7 AM to midnight. During working hours, at most of agents, but beginning in August, Internet (IN) customers
13regularagents,5Internetagents,and1shiftsupervisormay were served by a separate pool of agents. Thus, in August–
beworking. December,thecentercanbeconsideredtobetwoseparateser-
A simplified description of the path that each call follows vice systems, one for IN customers and another for all other
throughthe center is as follows. A customer calls one of sev- types. Second, as we discuss in Section 5, one aspect of the
eraltelephonenumbersassociatedwiththecallcenter,withthe service time data changed at the end of October. In several
number depending on the type of service sought. Except for instances, this article’s analyses are based on only the No-
rarebusysignals,thecustomeristhenconnectedtoaVRUand vember and December data. In other instances we have used
identifiesherself.WhileusingtheVRU,thecustomerreceives datafromAugust–December.Giventhechangesnotedearlier,
recordedinformation,bothgeneralandcustomized(e.g.,anac- thisensuresconsistencythroughoutthemanuscript.November
countbalance).It is also possible forthe customerto perform and December were also convenient, because they contained
someself-servicetransactionshere,and65%ofthebank’scus- noIsraeli holidays.In these analyses, we also restrictthe data
tomersactuallycompletetheirserviceviatheVRU.Theother to include only regular weekdays—Sunday–Thursday,7 AM–
35%indicatetheneedtospeakwithanagent.Ifanagentisfree midnight—becausethesearethehoursoffulloperationofthe
38 JournaloftheAmericanStatisticalAssociation,March2005
center. We have performed similar analyses for other parts of doesnotchangesignificantlywithinablock.Forconvenience,
thedata,andinmostrespectstheNovember–Decemberresults we used blocks of equal time length, L, although this equal-
do not differ noticeably from those based on data from other ity assumption could be relaxed. One can then consider the
monthsoftheyear. arrivals within a subset of blocks—for example, blocks at the
sametimeonvariousdaysorsuccessiveblocksonagivenday.
3. THEARRIVALPROCESS
Theformercasewould,forexample,testwhethertheprocessis
Figure 1 shows, as a function of time of day, the average homogeneouswithinblocksforcallsarrivingwithinthegiven
rate per hour at which calls come out of the VRU. These are timespan.
composite plots for weekday calls in November and Decem- Let T denote the jth ordered arrival time in the ith block,
ij
ber. The plots show calls according to the major call types. i=1,...,I.ThusT ≤···≤T ,whereJ(i)denotesthetotal
i1 iJ(i)
The volume of regular (PS) calls is much greater than that of numberofarrivalsintheithblock.ThendefineT =0and
i0
the other three types; hence those calls are shown on a sepa- (cid:1) (cid:1) (cid:2)(cid:2)
rate plot. [These plots were fit using the root–unroot method R =(J(i)+1−j) −log L−Tij , j=1,...,J(i).
described by Brown, Zhang, and Zhao (2001), along with the ij L−Ti,j−1
adaptivefreeknotsplinemethodologyofMaoandZhao(2003).
Undertheformalnullhypothesisthatthearrivalrateisconstant
Foramoreprecisestudyofthesearrivalrates,includingconfi-
within each given time interval, the {R } will be independent
dence and predictionintervals,see our Sec. 6 and also Brown ij
standardexponentialvariables,aswenowdiscuss.
etal.2001,2002a,b.]
LetU denotethejth(unordered)arrivaltimeintheithblock.
NotethebimodalpatternofPScall-arrivaltimesinFigure1. ij
It is especially interesting that IN calls do not show a similar Then the assumed constant Poisson arrival rate within this
bimodal pattern and in fact have a peak volume after 10 PM. block implies that, conditionally on J(i), the unordered ar-
(This peak can be partially explained by the fact that Internet rival times are independentand uniformly distributed, that is,
customers are sensitive to telephone rates, which significantly U i∼id U(0,L). Note that T = U . It follows that L−Tij
decrease in Israel after 10 PM, and that they also tend to be ariej independentbeta(J(i)+i1j −j,i1()j)variables[see,e.g.L,−pTrio,j−b1-
peoplewhostaylate.)
lem 6.14.33(iii) in Lehmann 1986]. A standard change of
variables then yields the conditional exponentiality of the R
3.1 ArrivalsAreInhomogeneous Poisson ij
given the value of J(i). [One may alternatively base the test
Common call center models and practice assume that the onthevariablesR∗ =j(−log Tij ),wherej=1,...,J(i) and
arrival process is Poisson with a rate that remains constant ij Ti,j+1
for blocks of time (e.g., half-hours),with a separate queueing Ti,J(i)+1=L. Underthe nullhypothesis,these will also be in-
dependentstandardexponentialvariables.]
modelfitted foreachblockoftime. A morenaturalmodelfor
The null hypothesis does not involve an assumption that
capturing changesin the arrivalrate is a time-inhomogeneous
the arrival rates of different intervals are equal or have any
Poisson process. Following commonpractice, we assume that
the arrival rate function can be well approximated as being otherprespecifiedrelationship.Anycustomarytest fortheex-
piecewiseconstant. ponential distribution can be applied to test the null hypothe-
Wenowconstructatestofthenullhypothesisthatarrivalsof sis.Forconvenience,weusethefamiliarKolmogorov–Smirnov
given types of calls form an inhomogeneous Poisson process test,eventhoughthismaynothavethegreatestpossiblepower
with piecewise constant rates. The first step in constructing againstthe alternativesof most interest. In addition,exponen-
our test involves breaking up the duration of a day into rela- tialQ–Qplotscanbeveryusefulinascertaininggoodnessoffit
tivelyshortblocksoftime,shortenoughsothatthearrivalrate totheexponentialdistribution.
(a) (b)
Figure1.ArrivalsinCalls/HourbyTimeofDay,WeekdaysinNovember–December.(a)PScalls;(b)IN,NW,andNEcalls.
Brownetal.:StatisticalAnalysisofaTelephoneCallCenter 39
Brownetal.(2002a)presentedquantileplotsforafewappli- Internetconsultationimposeson the system is greaterthan its
cationsof thistest. Forthe PSdata, we founditconvenientto shareintermsofpercentofcalls.Inearlierwork(Brownetal.
useL=6minutes.Fortheothertypes,weuseL=60minutes, 2002a)wealsoverifiedthatthefullcumulativedistributionsof
becausethesecallsinvolvedmuchlowerarrivalrates. theservicetimesarestochasticallyorderedinthesamefashion
We omit the plots here to save space and because they asthemeansinTable1.
demonstrateonlyminordeviationsfromtheidealstraight-line
4.1 VeryShortServiceTimes
pattern.OneexampleinvolvesarrivaltimesofthePScalls ar-
rivingbetween11:12AMand11:18AMonallweekdaysinNo- Figure2showshistogramsofthecombinedservicetimesfor
vember and December. A second example involves arrival of all types of service for January–October and for November–
IN calls on Monday, November 23, from 7 AM to midnight. December. These plots resemble those for PS calls alone,
Thiswasatypicalmidweekdayinourdataset. because the clear majority of calls are for PS. We see that in
Forbothoftheexamples,thenullhypothesisisnotrejected, thefirst10monthsoftheyear,thepercentageofcallswithser-
andweconcludethattheirdataareconsistentwiththeassump- vice<10secondswaslargerthanthepercentageattheendof
tion of an inhomogeneous Poisson process for the arrival of theyear(7%vs.2%).
calls.TherespectiveKolmogorov–Smirnovstatisticshaveval- Servicetimes<10secondsarequestionable.Andindeed,the
ues K = .0316 (p value≈.8 with n= 420) and K =.0423 manager of the call center discovered that short service times
(pvalue≈.9withn=172).Theseresultsaretypicalofthose were primarilycaused by agents who simply hungup on cus-
that we obtainedfrom various selections of blocksof the var- tomers to obtain extra rest time. (The phenomenon of agents
ious types of calls involving comparable sample sizes. Thus “abandoning” customers is not uncommon; it is often due to
overall, from tests of this nature, there is no evidence in this distorted incentive schemes, especially those that overempha-
datasettorejecta nullhypothesisthatthearrivalofcallsfrom size short average talk-time or, equivalently, the total number
theVRUisaninhomogeneousPoissonprocess. ofcallshandledbyanagent.)Theproblemwasidentified,and
As an attempt to further validate the inhomogeneousPois- stepsweretakentocorrectitinOctober1999.Forthisreason,
soncharacter,weappliedthismethodtothe48,193PScallsin inthelateranalysisofservicetimes,wefocusondatafromNo-
November and December in 6-minute blocks. With this large vemberandDecember.Suitableanalysescanbeconstructedfor
amount of data, one could expect to detect more than statis- theentireyearbyusingamixturemodelor,inasomewhatless
tically negligible departures from the null hypothesis because sophisticatedmanner,bydeletingfromtheservicetimeanaly-
of rounding of times in the data (to the nearest second) and sisallcallswithservicetimes<10seconds.
because arrival rates are not exactly constant within 6-minute
time spans. To compensate for the rounding,we “unrounded” 4.2 OnServiceTimesandQueueingTheory
the data before applying the test by adding independent uni-
Most applications of queueing theory to call centers as-
form (0,1) noise to each observation. (This unrounding did
sume exponentially distributed service times as their default.
noticeably improve the fit to the ideal pattern.) After the
The main reason for this is the lack of empirical evidence to
unrounding, the resulting Kolmogorov–Smirnov statistic was
the contrary, which leads one to favor convenience. Indeed,
K=.009.Thisisaverysmalldeviationfromtheideal;never-
modelswith exponentialservice times are amenableto analy-
theless,thepvalueforthisstatisticwithsuchalargen=48,963
sis,especiallywhencombinedwiththeassumptionthatarrival
isp≈.00007.[Toprovideanadditionalbenchmarkforevalu-
processesarehomogeneousPoissonprocesses.Thisistherea-
ating the (lack of) importance of this value, we note that this
son that M/M/N is the prevalent model used in call center
same statistic with n≈22,000would have had p value≈.05,
practice.
whichisjustacceptable.]
In more general queueing formulas, the service time often
affects performance measures through its squared coefficient
4. SERVICETIME
of variation, C2 =σ2/E2(S), E(S) is the average service time
s s
Thegoalofavisittothecallcenteristheserviceitself.Ta- and σ is its standard deviation. For example, a common use-
s
ble 1 summarizesthe mean,standarddeviation(SD),and me- ful approximationfor the average waiting time in an M/G/N
dianservicetimesforthefourtypesofserviceofmaininterest. model(Markovianarrivals,generallydistributedservicetimes,
Theveryfewcallswithservicetimes>1hourwerenotconsid- nservers),isgivenby
ered(i.e.,wetreatthemasoutliers).INcallshaslittleeffecton
(1+C2)
thenumbers.INcallshavethelongestservicetimes,withstock E[WaitforM/G/N]≈E[WaitforM/M/N]× s (1)
trading(NE)servicecallsnext.Potentialcustomers(NW)have 2
theshortestservicetime(whichisconsistentwiththenatureof (seeSze1984;Whitt1993).Notethatforlargecallcenters,this
thesecalls).Animportantimplicationisthattheworkloadthat formulamustbeusedwith care,asdiscussedbyMandelbaum
Table1.ServiceTimebyTypeofService,Truncatedat1Hour,November–December
Regularservice Potentialcustomers Internetconsulting Stocktrading
Overall (PS) (NW) (IN) (NE)
Mean 201 179 115 401 270
SD 248 189 146 473 303
Median 124 121 73 221 175
40 JournaloftheAmericanStatisticalAssociation,March2005
(a) (b)
Figure2.DistributionofServiceTime.(a)January–October(mean,185;SD,238);(b)November–December(mean,200;SD,249).
and Schwartz (2002). Thus average wait with general service densityasprovidedbyBrown andHwang(1993).Figure3(b)
timesismultipliedbyafactorof(1+C2)/2relativetothewait shows the lognormal Q–Q plot of service time. This does an
s
underexponentialservicetimes.For example,if servicetimes amazingly good imitation of a straight line. Nevertheless, the
are in fact exponential,then the factor is 1. Deterministic ser- Kolmogorov–Smirnovtestdecisivelyrejectsthenullhypothesis
vice times halve the average wait of exponential. In our data, ofexactlognormality.(TheKolmogorov–Smirnovstatistichere
theobservedfactoris(1+C2)/2=1.26. is K=.020.This is quitesmall, butstill muchlargerthan the
s
value of K=.009 that was attained for a similarly large sam-
4.3 ServiceTimesAreLognormal
plesizeintheinhomogeneousPoissontestofSec.4.)Weonly
Looking at Figure 2, we see that the distribution of service providethegraphstoqualitativelysupportourclaimoflognor-
times is clearly not exponential, as is assumed by standard mality.Thusthetruedistributionisveryclosetolognormal,but
queueing theory. In fact, after separating the calls with very isnotexactlylognormal.(Themostevidentdeviationisinthe
shortservicetimes,ouranalysisrevealsaremarkablefittothe lefttailofthehistogram,wherebothasmallexcessofobserva-
lognormaldistribution. tionsisevidentandtheeffectofroundingtothenearestsecond
Figure3(a)showsthehistogramoflog(servicetime)forNo- further interferes with a perfect fit.) This is a situation where
vemberandDecember,inwhichtheshortservicephenomenon averylargesamplesizeyieldsastatisticallysignificantresult,
wasabsentorminimal.Superimposedisthebestfittednormal eventhoughthereisno“practicalsignificance.”
(a) (b)
Figure3.Histogram(a)andQ–QPlot(b)oflog(servicetime),November–December.
Brownetal.:StatisticalAnalysisofaTelephoneCallCenter 41
After excludingshortservice times, the strongresemblance Y2i−1;T2i,Y2i}(cid:5)i=n/12(cid:6). The variance at T2i, σ2(T2i), is estimated
to a lognormal distribution also holds for all other months. byasquaredpseudoresidual,D2i,oftheform(Y2i−1−Y2i)2/2,
Italsoholdsforvarioustypesofcallers,eventhoughthepara- a so-called “difference-based” estimate. The difference-based
metersdependonthetypeofcall.Thismeansthatinthiscase, estimatorthatweusehereisasimpleonethatsufficesforour
amixtureoflognormalsisempiricallylognormal,eventhough purposes. In particular, our method yields suitable confidence
mathematically this cannot exactly hold. (See Mandelbaum intervalsforestimationofσ2. Moreefficientestimatorsmight
etal.2000,wherethephenomenonisdiscussedin thecontext improveour results slightly. There are many other difference-
of the exponentialdistribution.) Brown and Shen (2002) gave based estimators in the literature (see Müller and Stadtmüller
amoredetailedanalysisofservicetimes. 1987; Hall, Kay, and Titterington 1990; Dette, Munk, and
Lognormality of processing times has been occasionally
Wagner1998;Levins2002).
recognizedbyresearchersintelecommunicationsandpsychol- In the second step, we treat {T ,D }(cid:5)n/2(cid:6) as our ob-
ogy. Bolotin (1994) gave empirical results suggesting that 2i 2i i=1
served data points and apply local quadratic regression to ob-
the distribution of the logarithm of call duration is normal tain σˆ2(t ). Part of our justification is that under our model,
0
for individual telephone customers and a mixture of normals the {D }’s are (conditionally) independent given the {T }’s.
2i 2i
for “subscriber-line” groups. Ulrich and Miller (1993) and
A100(1−α)%confidenceintervalforσ2(t )isapproximately
0
Breukelen (1995) provided theoretical arguments for the log-
σˆ2(t )±z se (t ).
normality of reaction times using models from mathematical 0 α/2 σ2 0
Note that we use z , rather than a quantile from a chi-
psychology. Mandelbaum and Schwartz (2000) used simula- α/2
squareddistribution,asthecutoffvaluewhenderivingthefore-
tionstostudytheeffectoflognormallydistributedservicetimes
goingconfidenceinterval.Givenourlargedataset,thedegrees
onqueueingdelays.
offreedomarelarge,andachi-squareddistributioncanbewell
4.4 Regressionoflog(servicetimes)onTimeofDay approximatedbyanormaldistribution.
Theimportantimplicationoftheexcellentfittoalognormal 4.4.2 Estimation of ν(·). We now use µˆ(t0) and σˆ2(t0)
distributionisthatwecanapplystandardtechniquestoregress to estimate ν(t0), as eµˆ(t0)+σˆ2(t0)/2. Given that the estimation
log(servicetime)onvariouscovariates,suchastimeofday.For methods used for µ(t ) and σ2(t ), µˆ(t ) and σˆ2(t ), are as-
0 0 0 0
example, to model the mean service time across time of day, ymptoticallyindependent,wehave
we can first model the mean and variance of the log(service (cid:3) (cid:4) (cid:5)
time) acrosstime of day,then transformthe resultbackto the se µˆ(t0)+σˆ2(t0)/2 ≈ seµ(t0)2+seσ2(t0)2/4.
servicetimescale.[Shen(2002)gaveadetailedanalysisofser-
When the sample size is large, we can assume that µˆ(·)+
vicetimesagainstothercovariates,suchastheidentitiesofindi-
σˆ2(·)/2 has an approximately normal distribution. Then the
vidualagents(servers),aswellasreferencestootherliterature
corresponding100(1−α)%confidenceintervalforν(t )is
involvinglognormalvariates.] 0
(cid:5)
(cid:3)(cid:3) (cid:4) (cid:4)
Let S be a lognormally distributed random variable with
mean ν and variance τ2. Then Y =log(S) will be a normal exp µˆ(t0)+σˆ2(t0)/2 ±zα/2 seµ(t0)2+seσ2(t0)2/4 .
random variable with some mean µ and variance σ2. It is 4.4.3 Application and Model Diagnostics. In the analysis
wellknownthatν=eµ+σ2/2.Thisparameter(ratherthanµor thatfollows,weapplytheforegoingproceduretotheweekday
µ+σ2/2)istheprimitivequantitythatappearsincalculations callsinNovemberandDecember.Theresultsfortwointerest-
ofofferedload,asinSection7.Toprovideaconfidenceinterval ing service types are shown in Figure 4. There are 42,613PS
for ν, we need to deriveconfidenceintervalsfor µ and σ2 or, calls and 5,066 IN calls. To produce the figures, we use the
moreprecisely,forµ+σ2/2. tricubefunctionasthekernelandnearest-neighbortypeband-
For our call center data, let S be the service time of a call widths. The bandwidths are automatically chosen via cross-
andletT bethecorrespondingtimeofdayatwhichthecallbe- validation.
ginsservice.Let{S,T}n bearandomsampleofsizenfrom
i i i=1 Figure 4(a) shows the mean service time for PS calls as
thejointdistributionof{S,T}andsortedaccordingtoT.Then
i a function of time of day, with 95% confidence bands. Note
Y =log(S)willbethelog(servicetime)ofthecalls,whichare
i i theprominentbimodalpatternofmeanservicetimeacrossthe
(approximately)normallydistributed,conditionalonT.Wecan
i dayforPScalls.Theaccompanyingconfidencebandshowsthat
fitaregressionmodelofY onT asY =µ(T)+σ(T)(cid:4),where
i i i i i i thisbimodalpatternishighlysignificant.Thepatternresembles
(cid:4)|T areiidN(0,1).
i i thatforarrivalratesofPScalls(seeFig.1).Thisissuewasdis-
4.4.1 Estimation of µ(·) and σ2(·). If we assume that cussedfurtherbyBrownetal.(2002a).
µ(·) has a continuous third derivative, then we can use lo- Figure4(b)plotsananalogousconfidencebandforINcalls.
cal quadratic regression to derive an estimate for µ(·) (see OneinterestingobservationisthatINcallsdonotshowasim-
Loader1999).Supposethatµˆ(t ) isa localquadraticestimate ilar bimodal pattern. We do see some fluctuations during the
0
for µ(t ). Then an approximate 100(1−α)% confidence in- day,butthese are only mildlysignificant, giventhe wide con-
0
terval for µ(t ) is µˆ(t )± z se (t ), where se (t ) is the fidence band. Also notice that the entire confidence band for
0 0 α/2 µ 0 µ 0
standarderrorofthe estimate of the meanat t fromthe local INcallsliesabovethatofPScalls.Thisreflectsthestochastic
0
quadraticfit. dominancereferredtointhediscussionofTable1.
Our estimation of the variance function σ2(·) is a two- Standard diagnostics on the residuals reveal a qualitatively
step procedure. In the first step, we regroup the observa- very satisfactory fit to lognormality, comparable with that in
tions {Ti,Yi}ni=1 intoconsecutivenonoverlappingpairs {T2i−1, Figure3.
42 JournaloftheAmericanStatisticalAssociation,March2005
(a) (b)
Figure4.MeanServiceTime(PS)(a)versusTimeofDay[95%confidenceinterval(CI)],(b)MeanServiceTime(IN)versusTimeofDay
(95%CI).
5. WAITINGFORSERVICEORABANDONING covariates)can beestimated usingthe standardKaplan–Meier
product-limitestimator.
InSections3and4wecharacterizedtwoprimitivesofqueue-
One may plot the Kaplan–Meier estimates of the survival
ingmodels,thearrivalprocessandservicetimes.Ineachcase
functionsof R (timewillingto wait), V (virtualwaitingtime),
wewereabletodirectlyobserveandanalyzetheprimitiveun-
and W =min{V,R} (see Brownet al. 2002a).Thereis a clear
der investigation. We next address the last system primitive,
stochastic ordering between V and R in which customers are
customer patience and abandonmentbehavior, and the related
willing to wait (R) more than they need to wait (V). This
output of waiting time. Abandonment behavior and waiting
suggests that our customer populationconsists of patient cus-
timesaredeeplyintertwined.
tomers. Here we have implicitly, and only intuitively, defined
Thereisadistinctionbetweenthetimethatacustomerneeds
thenotionofapatientcustomer.(Tothebestofourknowledge,
to wait before reaching an agent and the time that a customer
systematicresearchonthissubjectislacking.)
is willing to wait before abandoning the system. The former
We also consider the survival functions of R for different
is referred to as virtual waiting time, because it amounts to
types of service. Again, a clear stochastic ordering emerges.
the time that a (virtual) customer, equipped with infinite pa-
Forexample,customersperformingstocktrading(typeNE)are
tience, would have waited until being served. We refer to the
latteraspatience.Bothmeasuresareobviouslyofgreatimpor- willingtowaitmorethancustomerscallingforregularservices
tance, but neither is directly observable, and hence both must (type PS). A possible empirical explanation for this ordering
beestimated. isthattypeNEneedstheservicemoreurgently.Thissuggests
A well-known queueing-theoretic result is that in heavily a practical distinction between tolerance for waiting and loy-
loadedsystems(inwhichessentiallyallcustomerswaitandno alty/persistency.
one abandons), waiting time should be exponentially distrib-
5.2 HazardRates
uted.(SeeKingman1962foranearlyresultandWhitt2002for
arecenttext.)Althoughoursystemisnotveryheavilyloaded, Palm(1953)wasthefirsttodescribeimpatienceintermsof
and in oursystem customersdoabandon,we findthat the ob- a hazard rate. He postulated that the hazard rate of the time
served distribution of time spent in the queue conforms very willingtowaitis proportionaltoa customer’sirritationdueto
well to this theoretical prediction (see Brown et al. 2002a for waiting. Aalen and Gjessing (2001) advocated dynamic inter-
furtherdetails). pretation of the hazard rate, but warned against the possibil-
itythatthepopulationhazardratemaynotrepresentindividual
5.1 SurvivalCurvesforVirtualWaitingTime
hazardrates.
andPatience
Wehavefounditusefultoconstructnonparametricestimates
Both times to abandonment and times to service are cen- of the hazard rate. It is feasible to do so because of the large
sored data. Let R denote the “patience” or “time willing to sample size of our data (about 48,000). Figure 5 shows such
wait” and let V denote the “virtual waiting time,” and equip plotsforRandV.
both with steady-state distributions. One actually samples The nonparametric procedure that we use to calculate and
W=min{R,V}, as well as the indicator 1{R<V}, for observ- plotthefiguresisasfollows.Foreachintervaloflengthδ,the
ing R or V. One considers all calls that reached an agent as estimateofthehazardrateiscalculatedas
censoredobservationsforestimatingthedistributionofR,and [numberofeventsduring(t,t+δ]]
vice versa for estimating the distribution of V. We make the .
[numberatriskatt]×δ
assumption that (as random variables) R and V are indepen-
dent given the covariates relevant to the individual customer. For smaller time values, t, the numbersatrisk and eventrates
Underthisassumption,thedistributionsofRandV (giventhe arelarge,andweletδ=1second.Forlargertimes,whenfewer
Brownetal.:StatisticalAnalysisofaTelephoneCallCenter 43
(a) (b)
Figure5.HazardRatesfor(a)theTimeWillingtoWaitforPSCallsand(b)VirtualWaitingTime,November–December.
areatrisk,weuselargerδ’s.Specifically,thelargerintervalsare time a customer is willing to wait to the mean time he or she
constructedtohaveanestimatedexpectednumberofeventsper needs to wait. The justification for calling this a “patience in-
intervalofatleastfour.Finally,thehazardrateforeachinterval dex”isthatforexperiencedcustomers,thetimethatoneneeds
isplottedattheinterval’smidpoint. towaitisinfactthattimethatoneexpectstowait.Althoughthis
Thecurvessuperimposedontheplottedpointsarefittedus- patienceindexmakessense intuitively,its calculationrequires
ing nonparametric regression. In practice we used LOCFIT the application of survival analysis techniques to call-by-call
(Loader1999),althoughothertechniques,suchaskernelproce- data.Suchdata maynotbe availableincertaincircumstances.
duresorsmoothingsplines,wouldyieldsimilarfits.Wechoose Therefore,wewishtofindanempiricalindexthatwillworkas
the smoothing bandwidth by generalized cross-validation. anauxiliarymeasureforthepatienceindex.
(We also smoothly transformed the x-axis, so that the obser- For the sake of discussion, we assume that V and R are in-
vationswouldbemorenearlyuniformlyplacedalongthataxis, dependentand exponentiallydistributed. As a consequenceof
beforeproducingafittedcurve.Wetheninverselytransformed theseassumptions,wecandemonstratethat
tthecehxn-iaqxuiessttohaittsvoarriigeidnathlefobramn.d)wWidethextopetraikmeenintetodawcictohufintttitnhge Patienceindex=(cid:7) mR = P(V<R).
m P(R<V)
V
increasedvarianceanddecreaseddensityoftheestimateswith
Furthermore,P(V<R)/P(R<V)canbeestimatedby(number
increasingtime.However,with ourdata, these techniqueshad
served)/(numberabandoned),andwedefine
littleeffect,andthuswedonotusehere.
Figure5(a)plotsthehazardratesofthetimewillingtowait (cid:7) numberserved
Empiricalindex= .
forPScalls.Notethatitshowstwomainpeaks.Thefirstpeak numberabandoned
occurs after only a few seconds. When customers enter the
Thenumbersofbothservedandabandonedcallsareveryeasy
queue, a “please wait” message, as described in Section 2, is
to obtainfromeither call-by-calldata ormore aggregatedcall
playedforthefirsttime.Atthispoint,somecustomerswhodo
center managementreports. We have thus derived an easy-to-
notwishtowaitprobablyrealizetheyareinaqueueandhang
calculate empirical measure from a probabilistic perspective.
up.Thesecondpeakoccursataboutt=60,aboutthetimethat
Thesamemeasurecanalsobederivedusingthemaximumlike-
the system plays the message again. Apparently, the message
lihoodestimatorsforthemeanofthe(right-censored)exponen-
increasescustomers’likelihoodof hangingupfor a brieftime
tialdistribution,appliedseparatelytoRandtoV.
thereafter, an effect that may be contrary to the message’s in-
Wecanuseourdatatovalidatetheempiricalindexasanes-
tendedpurpose(ormaybenot).
timate of the theoretical patience index. Recall, however, that
InFigure5(b),thehazardrateforthevirtualwaitingtimesis theKaplan–Meierestimateofthemeanisbiasedwhenthelast
estimated for all calls. (The picture for PS alone is very simi- observation is censored or when heavy censoring is present.
lar.)Theoverallplotrevealsratherconstantbehaviorandindi- Nevertheless, a well-known property of exponential distribu-
catesamoderatefittoanexponentialdistribution.(Thegradual tions is that their quantiles are just the mean multiplied by
general decrease in this hazard rate, from about .008 to .005, certainconstants,andweusequantileswhencalculatingthepa-
suggestsanissuethatmaymeritfurtherinvestigation.) tienceindex.Infact,becauseofheavycensoring,wesometimes
do not obtain an estimate for the median or higher quantiles.
5.3 PatienceIndex
Therefore,weusedfirstquartileswhencalculatingthetheoret-
Customer patience on the telephone is important, yet it has icalpatienceindex.
not been extensively studied. In the search for a better under- Theempiricalindexturnsouttobea verygoodestimate of
standingof patience,we havefoundarelative definitionto be thetheoreticalpatienceindex.Foreachof68quarterhoursbe-
ofuse.LetthemeansofV andRbemV andmR.Onecandefine tween 7 AM and midnight, we calculated the first quartiles of
the patience index as the ratio m /m , the ratio of the mean V andR fromthesurvivalcurveestimates.We thencompared
R V
44 JournaloftheAmericanStatisticalAssociation,March2005
¯
the ratio of the first quartiles to that of (number of served) Figure 1 presents a summary estimate of (cid:7)·(t), the average
to (number of abandoned).The resulting 68 sample pairs had of(cid:7)(t)overweekdaysinNovemberandDecember.
j
an R2 of .94 (see Brown et al. 2002a for a plot). This result For simplicity of presentation, here we treat together all
suggeststhatwecanusetheempiricalmeasureasanindexfor calls except the IN calls, because these were served in a sep-
humanpatience. aratesysteminAugust–December.Thearrivalpatternsforthe
Withthisinmind,weobtainthefollowingempiricalindices other typesof calls appearto be reasonablystable in August–
forregularweekdaysinNovemberandDecember:PS=5.34, December. Therefore, in this section we use the August–
NE = 8.71, NW = 1.61, and IN = 3.74. We thus find that December data to fit the arrival parameters. To avoid having
the NE customers are the most patient, perhaps because their to adjust for the short service time phenomenon noted in
businessis the mostimportantto them.On the other hand,by Section 4.1, we use only November and December data to fit
thismeasuretheINcustomersarelesspatientthanthePScus-
parametersforservicetimes.Also,hereweconsideronlyreg-
tomers. In this context, we emphasize that the patience index
ularweekdays(Sunday–Thursday)thatwerenotfullorpartial
measures time willing to wait normalized by time needed to
holidays.
wait. In our case (as previously noted), the IN customers are
Together, an arbitrary arrival rate (cid:7)(t) and mean service
inaseparatequeuefromthatofthePScustomers.TheINcus- time ν(t) at t define the “load” at that time, L(t)=(cid:7)(t)ν(t).
tomersonaveragearewillingtowaitslightlylongerthanthePS
Thisistheexpectedtimeunitsofworkarrivingperunitoftime,
customers(seeBrownetal.2002a).However,theyalsoneedto
aprimitivequantityinbuildingclassicalqueueingmodels,such
waitlonger,andoveralltheirpatienceindexislessthanthatof
asthosediscussedinSection7.
thePScustomers.
Briefly,supposethatoneadoptsthesimplestM/M/Nqueue-
Recall that the linear relationship between the two indices
ingmodel.Then,iftheloadisaconstant,L,overasufficiently
is established under the assumption that R and V are expo-
long period, the call center must be staffed, according to the
nentially distributed and independent. As Figure 5(a) shows,
model, with at least L agents; otherwise, the model predicts
however,thedistributionforRisclearlynotexponential.Sim-
thatthebacklogofcallswaitingtobeservedwillexplodeinan
ilarly, Figure 5(b) shows that V also displays some deviation
ever-increasingqueue. Typically, a manager will need to staff
fromexponentiality.Furthermore,sequentialsamplesofV are
notindependentofeachother.Thuswefindthatthelinearrela- the center at a√staffing level that is some function of L—for
example, L+c L for some constant c—to maintain satisfac-
tionissurprisinglystrong.
toryperformance(see Borst, Mandelbaum,andReiman 2004;
Finally,wenoteanotherpeculiarobservation:Thelinedoes
not have an interceptat 0 or a slope of 1, as suggested by the Garnettetal.2002).
foregoingtheory.Rather,theestimatedinterceptandslopeare
−1.82and1.35,whicharestatistically differentfrom0and1. 6.2 Independenceof(cid:7)(t)andν(t)
Weareworkingonprovidingatheoreticalexplanationthatac- In Section 5.4.4 we noted a qualitative similarity in the bi-
counts for these peculiar facts, as well as an explanation for modal pattern of arrival rates and mean service times. To try
the fact thatthe linear relationshipholdsso well, even though
to explain this similarity, we tested several potential explana-
the assumption of exponentiality does not hold for our data.
tions, including a causal dependence between arrival rate and
(The assumption of independence of R and V may also be
servicetimes.Wewereledtotheconclusionthatsuchacausal
questionable.)
dependenceis nota statistically plausibleexplanation.Rather,
weconcludedthattheperiodsofheaviervolumeinvolveadif-
6. PREDICTIONOFTHELOAD
ferentmixofcustomers,amixthatincludesahigherpopulation
This section reflects the view of the operations manager of
ofcustomerswhorequirelengthierservice.Thestatisticalevi-
a call center who plans and controls daily and hourly staffing
denceforthisconclusionisindirectandwasreportedbyBrown
levels. Prediction of the system “load” is a key ingredient in
et al. (2002a).Thuswe proceedunderthe assumption that ar-
thisplanning.Statistically,thispredictionisbasedonacombi-
rivalratesandmeanservicetimesareconditionallyindependent
nationoftheobservedarrivaltimestothesystem(asanalyzed
giventhetimeofday.
inSec.3)andservicetimesduringprevious,comparableperi-
ods(asanalyzedinSec.4).
6.3 CoefficientofVariation forthePredictionofL(t)
Inthediscussionthatfollowswedescribeaconvenientmodel
andacorrespondingmethodofanalysisthatcanbeusedtogen- Here we discuss the derivation of approximate confidence
erate predictionconfidenceboundsforthe load of the system. intervals for (cid:7)(t) and ν(t) based on observations of quarter-
More specifically, in Section 6.4 we present a model for pre- hourgroupingsofthedata.Theload,L(t),isaproductofthese
dicting the arrivalrate, and in Section 6.6 we presenta model two quantities. Hence exact confidence bounds are not read-
forpredictingmeanservicetime.InSection6.7wecombinethe ilyavailablefromindividualboundsforeachof(cid:7)(t)andν(t).
twopredictionstoobtainaprediction(withconfidencebounds) Asanadditionalcomplication,thedistributionsoftheindivid-
fortheloadaccordingtothemethoddiscussedinSection6.3. ual estimates of these quantities are not normally distributed.
Nevertheless, one can derive reasonable approximate confi-
6.1 DefinitionofLoad
denceboundsfromthecoefficientofvariation(CV)forthees-
In Section 3 we showed that arrivals follow an inhomoge- timateofL.
neousPoissonprocess.Welet(cid:7)(t)denotethetruearrivalrate For any nonnegative random variable W with finite posi-
j
of this process at time t on a day indexed by the subscript j. tivemeanandvariance,definetheCV(asusual)byCV(W)=
Brownetal.:StatisticalAnalysisofaTelephoneCallCenter 45
SD(W)/E(W). If U and V are two independentvariables and nearly precise even for rathe√r small values of λ. [One instead
W=UV,thenanelementarycalculationyields could use the si(cid:6)mpler form X or(cid:6)the version of Anscombe
(cid:5)
(1948)thathas X+ 3 inplaceof X+ 1;onlynumerically
CV(W)= CV2(U)+CV2(V)+CV2(U)·CV2(V). 8 4
smallchangeswouldresult.Ourchoiceisbasedonconsidera-
In our case, U and V correspond to (cid:7) and ν. Predictions tionsofBrownetal.(2001).]Additionally,V isasymptotically
for (cid:7) and ν are discussed in Sections 7.4 and 7.6. As noted normal(asλ→∞),anditmakessen(cid:6)setotreatitassuchinthe
earlier, these predictions can be assumed to be statistically modelsthatfollow.WethusletV = N + 1,andassumethe
jk jk 4
independent. Also, their CVs are quite small (<.1). Note model
that Lˆ(t) = (cid:7)ˆ(t)νˆ(t), and, using standard asymptotic nor- (cid:3) (cid:4)
m(cid:6)al theory, we can approximate CV(Lˆ)(t) as CV(Lˆ)(t) ≈ Vjk=θjk+εj∗k withεj∗ki∼idN 0,14 ,
CV2((cid:7)ˆ)(t)+CV2(νˆ)(t). θ =αβ +ε , (3)
jk j k jk
This leads to approximate 95% CIs of the form Lˆ(t) ±
2Lˆ(t)CV(Lˆ)(t). The constant 2 is based on a standard asymp- αj=µ+γVj−1,++Aj,
(cid:7)
toticnormalapproximationofroughly1.96. where Aj ∼ N(0,σA2), εjk ∼ N(0,σε2), Vj,+ = kVjk, and
6.4 Predictionof(cid:7)(t) Aj and εjk areindependentof eachotherandof valuesofVj(cid:9),k
(cid:9)
forj <j.Notethatα isarandomeffectinthismodel.Further-
j
BrownandZhao(2001)investigatedthepossibilityofmod- more, the model supposes a type of first-order autoregressive
elingtheparameter(cid:7)asadeterministicfunctionoftimeofday, structureon therandomdailyeffects.The correspondencebe-
dayofweek,andtypeofcustomer,andrejectedsuchamodel. tween (2) and (3) implies that this structure is consistentwith
Herewe constructa random-effectsmodelthatcanbeusedto anapproximateassumptionthat
predict(cid:7)andtoconstructconfidencebandsforthatprediction. (cid:1) (cid:9) (cid:2)
Themodelthatweconstructincludesanautoregressivefeature Rj= γ (cid:8) Nj−1,k+ 1 +Aj 2.
thatincorporatesthepreviousday’svolumeintotheprediction 4
k
oftoday’srate.
ThemodelisthusnotquiteanaturaloneintermsofR,butit
In the model, which we elaborate on later, we predict the j
appearsmorenaturalintermsoftheV in(3)andiscomputa-
arrivalonafuturedayusingarrivaldataforalldaysuptothat jk
tionallyconvenient.
day. Such predictions should be valid for future weekdays on
whichthearrivalbehaviorfollowsthesamepatternasthosefor Theparametersγ andβk needtobee(cid:7)stimated,asdoµ,σA2,
thatperiodofdata. and σε2. We impose the sid(cid:7)e condition βk2 =1, which cor-
Ourmethodofaccountingfordependenceontimeanddayis responds to the condition τk =(cid:6)1. The goal is then to de-
more conveniently implemented with balanced data, although riveconfidenceboundsforθjk= (cid:7)jk in(3),andsquaringthe
it can also be used with unbalanced data. For convenience, boundsyieldscorrespondingboundsfor(cid:7)jk.
we have thus used arrival data from only regular (nonholi- Theparametersin the model(3) caneasily be estimated by
day) weekdays in August–December on which there were no acombinationofleastsquaresandmethodofmoments.Begin
quarter-hourperiods missing and no obviousgross outliers in bytreatingthe{αj}’sasiftheywerefixedeffectsandusingleast
observedquarter-hourlyarrivalrates.Thisleaves101days.For squarestofitthemodel
eachday(indexedbyj=1,...,101),thenumberofarrivalsin V =αβ +(ε +ε∗).
each quarter hour from 7 AM–midnight was recorded as Njk, jk j k jk jk
k=1,...,68. As noted in Section 3, these are assumed to be This is an easily solved nonlinear least squares problem.
Poissonwithparameter(cid:7)=(cid:7)jk. It yields estimates αˆj, βˆk, and σˆ2, where the latter estimate is
One could build a fundamental model for the values of (cid:7) themeansquarederrorfromthisfit.Thenσ2 canbeestimated
ε
accordingtoamodeloftheform bymethodofmomentsas
Njk=Poisson((cid:7)jk), (cid:7)jk=Rjτk+εj(cid:9)k, (2) σˆ2=σˆ2− 1.
wheretheτ ’sarefixeddeterministicquarter-hourlyeffects,the ε 4
k
Rj’sarerandomdailyeffectswithasuitablestochasticcharac- Thenuse the estimates {αˆj} to constructthe least squaresesti-
(cid:9)
ter,andtheε ’sarerandomerrors.Notethatthismultiplicative matesoftheseparametersthatwouldbeappropriateforalinear
jk
structureisnatural,inthattheτk’splaytheroleoftheexpected modeloftheform
proportionoftheday’scallsthatfallinthekthinterval.Thisis
assumedtonotdependontheR’s,theexpectedoverallnumber αˆj=µ+γVj−1,++Aj. (4)
j
o(cid:7)fcallsperday.(Weaccordinglyimposethesideconditionthat Thisyieldsleast squaresestimates, µˆ and γˆ, andthe standard
τk=1.) meansquarederrorestimator,σˆ2,forthevarianceofA.
We instead proceed in a slightly different fashion that is A j
Theestimatescalculatedfromourdataforthequantitiesre-
nearlyequivalentto(2),butiscomputationallymoreconvenient
latedtotherandomeffectsare
and leads to a conceptually more familiar structure. The ba-
sisforourmethodisaversionoftheusualvariance-sta(cid:6)bilizing µˆ =97.88, γˆ =.6784
transformation.IfXisaPoisso√n(λ)variable,thenV= X+ 1 (withcorrespondingR2=.501), (5)
4
hasapproximatelymeanθ = λandvarianceσ2= 1.Thisis σˆ2=408.3, σˆ2=.1078 (becauseσˆ2=.3578).
4 A ε