Recycled incomplete identification procedures for blood screening

Recycled incomplete identification procedures for blood

screening

Citation for published version (APA):

Bar-Lev, S. K., Boxma, O. J., Kleiner, I., Perry, D., & Stadje, W. (2017). Recycled incomplete identification

procedures for blood screening. European Journal of Operational Research, 259(1), 330-343.

https://doi.org/10.1016/j.ejor.2016.10.005

Document license:

TAVERNE

DOI:

10.1016/j.ejor.2016.10.005

Document status and date:

Published: 16/05/2017

Document Version:

Publisher’s PDF, also known as Version of Record (includes final page, issue and volume numbers)

Please check the document version of this publication:

• A submitted manuscript is the version of the article upon submission and before peer-review. There can be

important differences between the submitted version and the official published version of record. People

interested in the research are advised to contact the author for the final version of the publication, or visit the

DOI to the publisher's website.

• The final author version and the galley proof are versions of the publication after peer review.

• The final published version features the final layout of the paper including the volume, issue and page

numbers.

Link to publication

General rights

Copyright and moral rights for the publications made accessible in the public portal are retained by the authors and/or other copyright owners and it is a condition of accessing publications that users recognise and abide by the legal requirements associated with these rights. • Users may download and print one copy of any publication from the public portal for the purpose of private study or research. • You may not further distribute the material or use it for any profit-making activity or commercial gain

• You may freely distribute the URL identifying the publication in the public portal.

If the publication is distributed under the terms of Article 25fa of the Dutch Copyright Act, indicated by the “Taverne” license above, please follow below link for the End User Agreement:

www.tue.nl/taverne Take down policy

If you believe that this document breaches copyright please contact us at: openaccess@tue.nl

providing details and we will investigate your claim.

ContentslistsavailableatScienceDirect

European

Journal

of

Operational

Research

journalhomepage:www.elsevier.com/locate/ejor

Innovative

Applications

of

O.R.

Recycled

incomplete

identification

procedures

for

blood

screening

Shaul K. Bar-Lev

a,∗

_{, Onno Boxma}

b

_{, Igor Kleiner}

a

_{, David Perry}

a

_{, Wolfgang Stadje}

c a Department of Statistics, University of Haifa, Haifa 31905, Israel

b EURANDOM and Department of Mathematics and Computer Science, Eindhoven University of Technology, P.O. Box 513, Eindhoven 5600 MB, The

Netherlands

c Institute of Mathematics, University of Osnabrück, Osnabrück 49069, Germany

a

r

t

i

c

l

e

i

n

f

o

Article history: Received 5 May 2015 Accepted 3 October 2016 Available online 7 October 2016 Keywords:

Group testing Blood screening Markov chain

Combinatorial urn problems Recycled group testing

a

b

s

t

r

a

c

t

Theoperationofbloodbanksaimsatthecost-efficientsupplyofuncontaminatedhumanblood.Eachunit ofdonatedbloodgoesthroughmultipletestingforthepresenceofvariouspathogenswhichareableto causetransfusion-transmitteddiseases.Thebloodscreeningprocessiscomprisedoftwophases.Atthe firstphase,bloodunitsarescreenedtogetherinpooledgroupsofacertainsizebytheELISA(Enzyme LinkedImmuno-SorbentAssay)testtodetectvariousvirus-specificantibodies.Thesecondphaseofthe screeningprocessisconductedbyPCR(PolymeraseChainReaction)testingoftheindividualbloodunits ofthegroupsfoundcleanbytheinitialELISAphase.

Thousandsofunitsofdonatedbloodarrivedailyatthecentralbloodbankforscreening.Each screen-ingschemehasassociatedtestingcostsand testingtimes.Inaddition,eachbloodunitarriveswithan expirationdate.Asaresult,theshorterthetestingtime,thelongertheresiduallifetimethatisleftfor thebloodunitforfutureuse.Thecontrollerfacesanaturalandwell-motivatedoperationsmanagement problem. Hewillattemptto shortenthetestingperiod and reducethetestingcostswithout compro-misingtoomuchonthereliability.Toachievethesegoals,weproposeanewtestingprocedurethatwe termRecycledIncompleteIdentificationProcedure(RIIP).InRIIP,groupsofpooledbloodunitswhichare foundcontaminatedintheELISAtestaredividedintosmallersubgroupsandagaingroup-testedbyELISA, andsoforth,untileventuallythePCRtestisconductedforthosesubgroupswhicharefoundclean.We analyzeandoptimizetheperformanceofRIIPbyderivingexplicitformulasforthecostcomponentsof interestandmaximizetheprofitassociatedwiththeprocedure.Ournumericalresultssuggestthatitcan indeedbeprofitabletodoseveralcyclesatELISA.

© 2016 Elsevier B.V. All rights reserved.

1. Introduction

The operation of blood banks, worldwide, is aimed at the supply of uncontaminated human blood. In the laboratories of theCentral Blood Services(CBS’s), each donated blood unit goes throughmultipletests. Theseareaimed atdetermining theunit’s blood type and the presence of various pathogens which are able to cause transfusion-transmitted diseases, such as Hepati-tis B (HBV), Hepatitis C (HCV), Human Immunodeficiency Virus (HIV) and Syphilis; see, e.g., Hourfar et al. (2008), Kantanen, Koskela, and Leinikki (1996), Kline, Brothers, Brookmeyer, Zeger andQuinn(1989),Monzonetal.(1992),Schottstedt,Tuma,Bünger, andLefèvre(1998),Steineretal.(2010),Strameretal.(2004).

∗ _{Corresponding author.}

E-mail addresses: barlev@stat.haifa.ac.il (S.K. Bar-Lev),o.j.boxma@tue.nl (O. Boxma), igkleiner@gmail.com (I. Kleiner), dperry@stat.haifa.ac.il (D. Perry), wstadje@uos.de (W. Stadje).

Until20yearsago,theroutinetestingwasdonewiththeELISA (Enzyme Linked Immuno-Sorbent Assay) test that detects virus-specific antibodies inthe blood. This test serves as a markerfor thenumberofantibodiesdetectedafterapersonisinfected. How-ever,itisimportanttotakeintoaccounttheeffectofthewindow period,definedastheperiodelapsingbetweenthe timeaperson is infectedby thevirus till a high concentration of antibodiesis developed.Thewindowperiodvariesfordifferenttypesofviruses. Examplesofaveragewindowperiodsforsomevirusesare:22days forHIV,70daysforHCVand60forHBV.Thismeansthatifa per-sonhasjustrecentlybeeninfectedby avirus,the ELISAtestfails todetectsuchaninfection.

Afewyearsago,anewtestcalledPCR(PolymeraseChain Reac-tion) has been developed. Ifa person wasjust recently infected, the PCR test immensely multiplies the number of antigens and thusmakesitpossibletodetectthepresenceofpathogensincases wherethe ELISA testfails todoso (for exampledueto the win-dowperiod effect).Hence the PCR testhas amuch higher sensi-tivity(probability ofa positive testresult foran infectedperson)

http://dx.doi.org/10.1016/j.ejor.2016.10.005 0377-2217/© 2016 Elsevier B.V. All rights reserved.

than theELISA test.However, PCR isvery expensivecomparedto ELISA.

Therefore, blood banks inthe USA, Israel andsome countries inEuropehaveestablishedascreeningprotocolwherebyallblood unitsare ELISA testedingroups andthosewhich testednegative forELISAarere-testedindividuallywithPCR(seee.g.Hourfaretal., 2008;Schottstedtetal.,1998;Strameretal.,2004).

The operation of blood bank systems is characterized by two crucial factors: (i) testing procedures and (ii) perishability. Test-ing is necessary, since only clean blood units are used forblood transfusion;groupsfoundcontaminated atELISA, andunitsfound contaminated at PCR, are discarded. Thousands of units of do-nated blood arrive daily at the central blood bank system for screening.Testing timesandtestingcostsare associatedwiththe screening process. In addition,each blood unit has an expiration date; after that it is perished. The controller faces a natural and well-motivatedoperationsmanagementproblem.Ontheonehand, thereistheneedtomake thetestingperiodasshortaspossible; ontheotherhand,carefultestingisrequired,whichtakestimeand iscostly. Thisraises theneedto findmoreefficientgroup testing procedures,withtherestrictionofincompleteidentification.

Since ELISA testing is relatively cheap, we propose a new screeningprocessthatwetermRecycledIncompleteIdentification group testing Procedure (RIIP), by which groups ofpooled blood units, which are found contaminated in the previous ELISA cy-cle, aredividedintosmallersubgroups andagaingroup-testedby ELISA,andsoforth,untilfinallythePCRtestisconductedforthose subgroupswhicharefoundclean.

The goal of this paper is to provide an analysis and opti-mization oftheperformance ofRIIP, inparticularminimizing the costs (or maximizing the profit) associated with the test proce-dure.Somewhatrelatedpapersare(Bish,Bish,Xie,&Slonim,2011; Xie,Bish, Bish,Slonim, &Stramer,2012). Thosepapers alsostudy theproblemofselectinganeffectivesetofscreeningtestsfor do-natedblood.Theyfocusontheproblemofminimizing transfusion-transmitted infectious diseases, under certain budget constraints (Bish et al., 2011) and in addition also under waste constraints (Xie et al., 2012). Other optimization problems regarding blood management that havereceived much attentionare supplychain management and inventorymanagement. We refer to Beliën and Forcé (2012)for avery thoroughliterature reviewof inventoryand supply chainmanagement ofblood products, and tothe two re-cent studies (Blake & Hardy, 2014) andCivelek, Karaesmen, and Scheller-Wolf(2015)forinterestingnewresultsoninventory man-agementpolicies.

In thenext subsectionswe describe some features ofa blood banksystem. InSection 1.1wedescribe theseparationprocessof each blood unit into its three components along with their ex-piration dates andassociated costs. The concept ofgroup testing is reviewed in Section 1.2. In Section 1.3 we present the group testingprocedures(complete andincompleteidentification proce-dures)thatare currentlyinusebybloodbanks.Thefurther orga-nizationofthepaperisoutlinedinSection1.4.

1.1. Bloodcomponents

Blood consists of several components: Red blood cells (RBC), plasmaandplatelets.Processingthewholebloodunitsintothe dif-ferentcomponents isdone inparallel tothetestingstage. Whole blood units are separated intodifferent components, whichhave different biological functions, storage conditions and expiration dates.Theywillbesuppliedtodifferentpatientsaccordingtotheir medicalneeds:

1. RBC – This component, which is separated from the whole blood unit within 8–24 hours from collection, can be used

within 35–42 days, depending on whether additive solution is added. Mostly, the cost of an unmodified packed RBC unit for the hospitalsis $40 andthat of a leukodepleted(solution added)unitis$70.

2.Plasma – Plasma units are automatically made upon the pro-duction of an RBC unit. The corresponding cost for a plasma unitisaround$40;hospitalsacquire28%ofalltheplasmaunits produced.

3.Platelets– Fromeachwholebloodunitonerandomplateletunit is separated,which can be used for atmost5 days. The cost of producingrandom plateletunitsis about$40/unit(Bar-Lev, Stadje,&VanderDuynSchouten,2005).

Withrespecttotheaboveexpirationdates,thefollowingshould betakenintoaccount.Onaverageittakesabout15hourstillblood samplesarrive in the CBSafter the momentthey havebeen do-nated. The average processing time of an ELISA test is around 1 hourandcosts around$2.5 peraverage group(of bloodunits) of size10,whereas thatofthePCR isaround6hours withan asso-ciatedcostofabout$85per bloodunit.Suchtime constraintsare vitalforplatelets’shelf-life,butlesssignificantforRBC.However, such processingtimesshouldbe takeninto accountin anyblood screeningprocedure.

1.2.Grouptesting

Theissueofblood transfusionmightbe a questionoflife and death.Thismeansthatit isofparamountimportance thatall the blood unitsthat enter the shelf are clean.Therefore, a necessary requirementisameticulousinspectionofallthebloodunits. How-ever,sincethousandsofbloodunits(bloodsamples) arriveatthe centralbloodbankevery day,anaturalscreeningproceduremust be based on the idea of group testing; otherwise, the screening processwilltaketoolongandthecostsofthisprocesswillbetoo high.

Grouptestingdealswiththeclassificationoftheitemsofsome populationintotwocategories:‘good’and‘defective’.Itisassumed thattheitemsaregrouptestable,i.e.,foranysubset ofthe popu-lationit is possibleto carry out a simultaneoustest (grouptest) withtwopossibleoutcomes:‘success’(alsocalled‘clean’,or ‘nega-tive’),indicatingthatallitemsinthesubsetaregood,and‘failure’ (also called ‘contaminated’, or ‘positive’), indicating that at least one of the items in the subset is defective – without knowing whichor how manyare defective. Acontaminated group can be subject to further screening, or be scrapped. Employing suitably designedprocedures ofthis kind leadsto a significant reduction ofthenumberofrequiredtestsandthusofscreeningcosts,under controlled probabilities of misclassification. A group testing pro-cedure isthereforea cost-efficienttechnique. It hasbeenapplied invarious areas,firstandforemost, forblood screeningto detect various viruses, forDNA screening,as well as forquality control forindustrial productionsystems (see,e.g.,Bar-Lev,Blanc,Boxma, Janssen,&Perry,2013;Bar-Lev,Boneh,&Perry,1990;Bar-Lev, Par-lar,&Perry,1995;Bar-Lev,Parlar,Perry,&VanderDuynSchouten, 2007; Bar-Lev, Stadje, & Van der Duyn Schouten, 2003; 2004; 2005;2006;Beliën&Forcé,2012;Bishetal.,2011;Blake&Hardy, 2014;Chick,1996;Civeleketal.,2015;Du&Hwang,2000;Feller, 1968; Gastwirth & Johnson, 1994; Hammick & Gastwirth, 1994; Hanson,Johnson,&Gastwirth,2006;Hourfaretal.,2008; Hughes-Oliver & Rosenberger, 2000; Johnson & Gastwirth, 2000; Litvak, Tu, & Pagano, 1994; Macula,1999a; Macula, 1999b;Tu, Litvak,& Pagano, 1995; Uhl, Liu, Walther, Hess, & Naiman, 2001; Wein & Zenios,1996;Wolf,1985;Xie,Tatsuoka,Sacks,&Young,2001;Xie etal., 2012; Yamamura & Ishimoto, 2009; Zhu, Hughes-Oliver, & Young,2001).

1.3.Completeversusincompleteidentification

One maycurrentlydistinguishtwo typesofidentification test-ingprocedures:complete andincomplete.Thepurpose ofa com-plete identificationgroup testing procedureis to classify each item ina givenpopulationaseithercleanordefective.Thisisdoneby testinggroupsofsize m (adecision parameter) inthe ELISA sta-tiononly.Ifagroupisfoundcleanitisaggregatedforblood trans-fusionpurposes,otherwise,ifitisfoundcontaminateditisfurther re-testedbydividingitintosubgroups.Suchaprocedurecontinues tilleachiteminagivenpopulationisappropriatelyclassified.

One could imagine two managerialreasons inwhich complete identification procedures are inefficient. The first one is that the testsaretoo expensive. Then acomplete identificationprocedure leadstoahighexpectednumberoftestsand, asaresult, tohigh expectedtotaltestingcosts.Thesecondoneisthattheshelf-lifeis short(recallthattheshelf-lifeofplateletsisatmost5days).Then thehigherthenumber oftests,the shortertheresidual shelf-life ofcleanitemsontheshelf.

The idea of an Incomplete Identification group testing Proce-dure(IIP) was first introduced in Bar-Lev et al. (1990) for some industrial problem and was subsequently further developed for bloodscreening(e.g.Bar-levetal., 2013;Bar-Levetal.,1995; Bar-Levet al., 2007; Bar-Lev et al., 2003, 2004, 2005, 2006). An IIP startsasabovebyfirsttestinggroupsofsizemintheELISAstation. However,asopposedtothepreviouscase,agroupfound contam-inatedisscrapped;otherwise,itisaggregatedandsenttothePCR stationforindividualtesting.Suchaprocedureiscost-wiserather efficientasitsignificantlydecreasesthenumberoftests(whether groupedorindividual).It isparticularlyefficientwhenthe preva-lencerateofthe“deficiency” (like theprevalencerateof,say,HIV inthepopulation)israthersmall.

In order to give some idea how real data are processed we mentionthefollowing.InWesternEuropean countries(aswell as inthe US,Israelandsome other countries)about35,000–50,000 blooddonationsareneededper1millionpersonsperyear.The fol-lowingdatafortheyears2011–2012havebeenprovidedtousby theIsraeliCentralBloodBank.Thedatadescribeperyearthe num-ber of blood donations(blood units), the number of blood units foundcontaminated attheELISA station(includingHIV, HBVand HCV) and the number of units found clean at the ELISA station butthenfoundcontaminatedatthePCRstation(forthetwoyears 2011and2012 thepopulation size ofIsrael wasabout7,800,000 and7,900,000,respectively).

Donations Year Cases confirmed positive by ELISA HBV HCV HIV Total 294,117 2011 126 62 13 201 298,470 2012 143 62 3 208

The table suggeststhat the estimatedprobability offinding a givenunittobecontaminatedattheELISAstationisapproximately 0.00068 (lateron such a probability will be denoted by

ε

). Note thatthisestimatedprobabilityissmallerthantheprevalencerate ofthecontaminatingvirus(es)inthepopulationasthoseinfected persons who are aware of their situation usually do not donate bloodsamples.Thismeansthatiftheprobability ofcontaminated blood units is not known, for some reason, then the prevalence rateinthe populationcanbe usedasan upperbound forsucha probability.Allofthe positiveELISA unitswere alsojudged posi-tivebyPCR, but,inaddition,therewere another12unitsin2011 and 13 in 2012 that were only found positive by PCR (but not by ELISA). Hence, the estimated probability of those units found cleanbyELISAbutthenfoundcontaminatedbyPCRis0.00004for both years 2011and 2012 (laterthis probability will be denoted by

γ

).

1.4. Organizationofthepaper

Thepaperisorganizedasfollows.InSection2wedescribethe RIIPinmoredetail,includingcosts andtimesinvolvedaswellas related stopping times. We also define an appropriate objective function, aiming at maximizing the profit of the RIIP operation.

Sections3and4aredevotedtothederivationofthedistributional behavioroftherelevantstoppingtimesandallthefunctionals oc-curringintheobjectivefunction.Weintroduce anunderlying un-observablenonhomogeneousMarkovchain,whosedistributionwill be determined in closedform, and show that all distributions of interestandthusthecompleteobjectivefunctioncanbeexpressed in terms of this Markov chain. It will turn out that the result-ingformulasforthesedistributionsarerathercomplicated,making theobjectivefunctionquiteintricate andatpresentimpossibleto handlenumericallywithordinarycomputers.Thereforewedevote

Section 5to introducingapproximationsoftherelevant function-als.Theseapproximationseasethenumericalevaluationofthe ob-jective function. Numerical examples are presented in Section 6, along with some sensitivity analysis of the involved parameters anddecisionvariables.Section7containsconclusions.

2. TheRIIPmodel:description,stoppingtimes,parametersand objectivefunction

In the RIIP, blood samples are tested at two consecutive sta-tions:firsttheELISA stationandsubsequentlythePCR station.In thesequel,weconsidertestingonaweeklybasis.Aninitialsupply ofnbloodunitsperweek(accordingto thedatainSection 1.3,n

wouldbe around 6000)is divided intok1 groupsofsize m1; m1

isadecisionvariable.Thesegroupsareinitiallygroup-screenedat theELISAstation.IfagroupisfoundcleanitissenttothePCR sta-tionforindividualscreening(i.e.,eachofthem1unitsinthegroup

isscreenedindividually).Otherwise,ifthegroupisfound contam-inated, the m1 units inthe contaminated group are divided into

k2 subgroups ofsize m2,all ofwhich are recycledand resent to

the ELISA station for furthergroup-screening. The whole process then repeats itself for the next recycle (each contaminated sub-groupisdivided intok3 subgroupsofsizem3,whicharerecycled

andresent to the ELISA station forfurther group-screening), and so forth.The reasonforrecycling attheELISA station isthe cost reduction, asthe testing cost atthe ELISA station is significantly smaller($2.5)thanthatatthePCRstation($85per bloodunit, cf.

Section1.1).SucharecyclingattheELISAstationalsomakessense froma time perspective,as theprocessingtime atthis stationis around1hourpergroupwhereasthatatthePCRstationisaround 5–6hours.The aimofthetestingistosatisfy agivendemandof

dbloodunitsinacost-efficientmanner.Forthisonehastodecide onthe numberofcyclesattheELISA stationandthegroup sizes

mi.It shouldbenotedthatrepeatedgrouptestingdoesnotaffect

thequalityofthescreenedblood.

Inthesequel,weshallsaythatthebloodtestingprocesshasr cyclesifthelastsubdivisionisinkrgroupsofsizemr – soifone

doesonlyonegrouptest,thenr=1.Theefficientnumberofcycles dependsonseveralfactors.

First,thegroup sizessatisfy m1 > m2 >..., andaftera certain

numberofcycles,thenumberofitemsintheresultingsubgroups becomes too small, making further recycling no longer worth-while.Accordingly,therecyclingisstoppedatthesmallestmi≥ ˆm,

wheremˆ is the smallestpermissible group size.There are also a lowerboundmˆ andanupperboundm˜ onthenumberofsamples that can be pooledtogether dueto technicalrestrictions, so that

ˆ

m≤m1≤ ˜m;inpracticeusuallynotlessthan4andnomorethan

64samples aretakenin onegroup. Second, theprocessingtimes causean upper bound onthe number ofcycles. Werequire that thetotalprocessingtimeofall cyclesattheELISAstationwillnot

Fig. 1. Flowchart of the group testing procedure.

exceedapredeterminedtimet0.Weassumethat everytest takes

a time interval offixed (constant) length,regardless ofthe batch size; every test takesTelisa time units. Then the testingis

termi-natedatthelatestafterthe(t0/Telisa)thtest.

Third, we assume that there isa prespecified upperbound c0

forthetotal costoftheELISA cyclesandthat thecostper testis constant,sayCelisa.Therefore,thecostlimitisreachedafterc0/Celisa

tests.

Notethatasthelargestgroupsizedoesnotexceed64whereas thesmallestoneis4,itfollowsthatthelargestvalue ofris5

(

= log₂64_{− 1}

)

cycles.

LetNh bethetotalnumberoftestsafterhcycles.Accordingto

the above constraints,the ELISA recyclingprocess stops afterthe

τ

thcyclewhere

τ

=. sup



h:mh≥ ˆmandNh≤ t0 Telisa andNh≤ c0 Celisa



. (1)

Wecanequivalentlywritethisstoppingtimeas

τ

=min[

τ1

,

τ

C], (2) where

τ1

=. max

{

h:mh≥ ˆm

}

, (3) and

τ

C=. max

{

h:Nh≤ C

}

, (4) with C:=min



t0 Telisa , c0 Celisa



. (5)

Fig.1displaysaflowchartoftheproposedgrouptesting proce-duredescribedinthissection.

Wenowformulateaprofitobjectivefunction,whichhastobe maximizedwithrespecttothedecisionparametersm1>m2>...

Thefollowingnotationwillbeusedthroughout.

• Probabilityparameters

ε

-theprobabilitythataunitisfoundcontaminatedbyELISA.

γ

- the conditional probability that a unit is found contami-natedbyPCRgivenitwasfoundcleanbyELISA.

• Countingrandomvariables

Thevariablesareconsideredonaweeklybasis.

M-thetotalnumberofgrouptestsattheELISAstation.

N-thetotalnumberofunits(outofn)foundcleanbyELISA.

N∗ -the totalnumberofunits (outofn) foundclean by both ELISAandPCR.NotethatN∗≤ N≤ nand

E

(

N∗

)

=

(

1−

γ

)

E

(

N

)

. (6)

• Costs,penaltiesandrewards

Cpur -the purchasingcost for asingle unit (relatedto

acquir-ingtheninitialunits).Thisisthecostcausedbycollecting blood donationsfrom various locations(and thus is composed ofcosts

oftechnicians,nurses,vehicles,transportationandthelike)and/or paymentssometimesgiventodonors(asisthecaseintheUS).

Celisa -the costofa testfora batch ofarbitrarysizem atthe

ELISAstation.Accordingtopractitioners(asProfessorEilatShinar, anhematologist,thedirector oftheIsraeli CBS),theELISA cost is hardlyaffectedbythebatchsizeaslongasitdoesnotexceed64.

Cpcr-thecostoftestingasingleunitatthePCRstation.

Cpenalty-thepenaltycostfornotsatisfyingademandfora

sin-gleunit whichhasbeentestedclean.Penaltycosts occurwhena hospitalwhosedemandisnotsatisfiedusesitsownsources(e.g.a localbloodbank)tosatisfytherequireddemand.

Rw -therewardforasatisfieddemandunit.

ˆ

Rw -therewardforany‘surplus’unit beyondtherequired

de-mand.

NotethatRw istheprice(perbloodunit)acentralbloodbank

is paid by hospitals. In case of ‘surplus’, the units are kept but mightperishifnodemandoccurs.HenceRw>Rˆw.

Thisyieldsthefollowingcomponentsoftheobjectivefunction:

Costcomponents

Cpurn-thetotalcostforpurchasingnunits.

CelisaM-thetotalcostfortestingMgroupsattheELISAstation.

CpcrN-thetotalcostfortestingNunitsatthePCRstation.

C_penalty[d_{− N}∗]I_(N∗_≤d) - the total penalty for not satisfyingthe full demand d(here andin the sequel, I(· ) denotes an indicator

function).

Rewardcomponents

RwN∗I₍N∗≤d)-thetotalrewardforthesatisfied(partofthe)

de-mand. ˆ

Rw

(

N∗− d

)

I(N∗_>d)-thetotalrewardforthedemandsurplus. Thusthetotalrewardisgivenby

RwN∗I(N∗_≤d)+



Rwd+Rˆw

(

N∗− d

)



I₍N∗_>d). (7)

Combining all of the above, the (random) profit P=

P

(

m1,m2,...

)

associatedwiththeprocedureis

P=RwN∗I₍N∗_≤d)+



Rwd+Rˆw

(

N∗− d

)



I₍N∗_>d)

−



Cpurn+CelisaM+CpcrN+Cpenalty[d− N∗]I(N∗_≤d)



. (8)

The objectivefunctionisthen givenby theexpectedprofitP˜₌ E[P

(

m1,m2,...

)

]: ˜ P=E[RwN∗I₍N∗_≤d)]+E



{

Rwd+Rˆw

(

N∗− d

)

}

I₍N∗_>d)



−



Cpurn+CelisaE[M]+CpcrE[N]+CpenaltyE[

(

d− N∗

)

I(N∗_≤d)]

. (9) Werewritethisas ˜ P=P1+P2− P3− P4− P5− P6, (10) where

⎧

⎪

⎪

⎪

⎪

⎨

⎪

⎪

⎪

⎪

⎩

P1=E[RwN∗I₍N∗_≤d)], P2=E



{

Rwd+Rˆw

(

N∗− d

)

}

I(N∗_>d)



, P3=Cpurn, P4=CelisaE[M], P5=CpcrE[N], P6=CpenaltyE[

(

d− N∗

)

I(N∗_≤d)]. (11)

Tomaximize thisobjective function we haveto find explicit ex-pressions for all its ingredients, i.e., the expected values on the rightside of (9). Thisis carried out in the next two sectionsby derivingtheunderlyingdistributionsinclosedform.

3. TheunderlyingMarkovchain

InthissectionweintroduceaMarkovchain(MC)thatcaptures theessenceofthe testingprocedure.The evolutionof thisMC is basedonthefullinformationonthenumberofcontaminatedunits

intheformedgroupsineachstageofthetestingprocess,sothatit isonlypartiallyobservable.Wedetermineitsdistributioninclosed form, and we show that thedistributions ofall quantities of in-terest (M,N,N∗ featuring inthe profitobjectivefunction)can be expressedintermsofthisMC. ThereforetheMCserves asa con-venienttooltodeterminetheobjectivefunctioninclosedform.

Recallthat,afterconditioningontheinitialsupplyofthe num-berofcleanitemsinthepopulationofsizen,thebasicRIIPmodel can be describedas follows. Among n itemsr are contaminated andtherestisclean.In thefirststage theitemsaresplit intok1

groups,eachofgroupsizem1=n/k1,whicharepooledandtested

together.Eachofthegroupsfoundcontaminatedisthensplitinto

k2 subgroups,each ofsizem2=m1/k2, whicharethentested

to-gether.Theprocedure isiterated.Letusfirstassume that itis it-erated until in thefinal stage only contaminated ‘groups’ of size 1 are left, which are then tested and the complete picture be-comesknown.Laterwewilltruncatethisprocedureanduse stop-ping timesoftheform(2),butwe willseethatthismodification caneasilybetakencareof.Letpbe thenumberofcycles,sothat

n₌p_j=1k_j,and

m1>· · · >mp−1>mp>1.

Notethatk1=n/m1andkj=mj−1/mj for j=1,...,p (in

particu-lar,thek_jhavetobechosensuch thatk1 dividesnandkjdivides

mj−1forj>1).

We assume that the events {ith item is clean}, i=1,2,...,n,

areindependent.Thisassumptionisreasonableiftheitemscome fromdifferentdonors.Iftheinitialsupplycontainsitemsdonated bythesameindividual,onlyoneitemperindividualistakentobe tested. Itfollowsthat inourmodelthenumberofcleanitemsin thepopulationisbinomiallydistributed.

For j=1,2,...,p letY_ljbethenumberofgroupsinthejth cy-clethat containexactly lcontaminated items. For j=0let Y0 _be

the initialnumber ofcontaminated items. We introducethe ran-domvectors

Yj=

(

Y₀j,Y₁j,...,Ymjj

)

, j=0,1,2,...,p

wherem0=0andthusY00=Y0.ItisimportanttonotethattheY

j l,

l≥ 1,arenotobservableduringthetestingprocess.NotethatY₀jis thenumberofgroupsfoundcleaninthejthcycleandmj

l=1Y

j l is

thenumberofgroupsfoundcontaminatedinthejthcycle.OnlyY₀j

andm_l=1j Y_ljbecomeknownafterthejthcycle,butnotY₁j,...,Ymj_j

individually.

Wewill seethatall distributions ofinterestto usandall cost andrewardfunctionalswemaywanttoconsidercanbeexpressed intermsofthejointdistributionof

(

Y0_,Y1_,...,Yp

₎

_,

which isa sequence of random vectorsof successivedimensions 1,m1+1,...,mp+1. This finite sequence is a nonhomogeneous

Markovchain,sotodetermineitsjointdistributionweonlyhaveto deriveallitsone-steptransitionprobabilities.(Notethatthe distri-butionofY0 _{isbinomial.)Thesetransitionprobabilitiescanbe}

de-rivedbymeansofsomecombinatorialconsiderationstowhichwe turnnow.

3.1. Acombinatorialurnproblem

We need to solve the following auxiliary urn problem, which apparently has not been treatedin the voluminous literature on urnmodels(seee.g.Kolchin,1978;Kotz&Johnson,1977).Consider an urn containing rred ballsand n− r white balls. Take out all ballsinkgroupsofequalsizem=n/kaccordingtothe equidistri-bution.Suchequidistributedgroupscanforexamplebegenerated asfollows:

(1) Numbertheballsby1,...,n.

(2) Take a uniformly distributed random permutation

(τ

1,...,

τ

n

)

onthesetofallpermutationsof

{

1,...,n

}

.

(3) Take as groups the sets

{

τ

1,...,

τ

m

}

,

{

τ

m+1,...,

τ

2m

}

,...,

{

τ

(k−1)m,...,

τ

n

}

.

Forl=0,1,2,...letYl be thenumberofgroupscontaining

ex-actlylredballs;note thatYl=0forl>r.Theproblemisto

com-putetheprobability

p

(

y0,...,yr |n,r,m

)

=P

(

Y0=y0,...,Yr=yr

)

.

Solution:Weonlyconsiderthosevaluesofy0,...,yrforwhich

the probability is positive. Thus the yi’s are nonnegative integers

satisfying

y0+· · · +yr=k, (12)

y1+2y2+3y3+· · · +ryr=r, (13)

my0+

(

m− 1

)

y1+

(

m− 2

)

y2+· · · +

(

m− r

)

yr=n− r. (14)

Theseequationsexpressthattherearek groupsintotal,thetotal numberofredballsisrandthatofwhiteballsisn_{− r.}

Ofcourse p

(

y0,...,yr

|

n,r,m

)

isproportionaltothenumberof

ways,sayh

(

y0,. . .,yr

|

n,r,m

)

,tosplitthenballsaccordingly.

Inafirststep,toachieveY0=y0,wehavetotakeouty0groups

ofsizemfromthen− r whiteballs.Forthisthereare h0

(

y0 |n,r,m

)

= 1 y0!



n−r m



n−r−m m



· · ·



n−r−

(

y0−1

)

m m



(15)

possibilities;the factor1/y0! occurssincewe want an unordered

setofgroups.

Howmanypossibilitiesarethere toachieveYl=yl,giventhat

Y0=y0,...,Yl−1=yl−1? Wehavetoformylgroupsofsizemeach

containinglredballsandm− l whiteballs,wherethewhiteballs haveto bechosen fromthen− r− y0m− · · · − yl−1

(

m− l+1

)

re-maining whiteballsandtheredballshavetobechosen fromthe

r− y1− 2y2− · · · −

(

l− 1

)

yl−1 remaining red balls. The unordered

numberofwaysforthisis hl

(

yl |y0,...,yl−1,n,r,m

)

= 1 yl!



r− y1− 2y2− · · · −

(

l− 1

)

yl−1 l



×



n− r− y0m− ...− yl−1

(

m− l+1

)

m− l



× · · · ×



r− y1− 2y2− · · · −

(

l− 1

)

yl−1− l

(

yl− 1

)

l



×



n− r− y0m− · · · − yl−1

(

m− l+1

)

−

(

yl− 1

)(

m− l

)

m− l



. (16) Therefore, h

(

y0,. . .,yr |n,r,m

)

=h0

(

y0 |n,r,m

)

r  l=1 hl

(

yl |y0,...,yl−1,n,r,m

)

.

Notethatthetermsinformula(16)aresymmetricwithrespectto the red and white balls(they also consider all possible arrange-mentsofwhiteballswithinagroupandallpossiblearrangements of red balls). This also holds for (15), as becomes obvious after multiplyingtheright-handsideof(15)by

(



r_r



)

y0_.

Formula(16) greatlysimplifies after cancelingout several fac-tors.Notice thatthefirsttermof(16)(fory1)startswitha factor

(

n− r− y0m

)

!, which cancels against the

(

n− r− y0m

)

! term in (15);noticethatonecanrewrite(15)as

h0

(

y0

|

n,r,m

)

=

1

y0!

(

n− r

)

!

(

m!

)

y0

(

n− r− y₀m

)

!.

Inthisway,a lotoftermscancel whenwe considerh₌r_j=0h_j,

andweget: h

(

y0,...,yr

|

n,r,m

)

=

(

n− r

)

!r! y0!...yr! 1

(

m!0!

)

y0

((

m− 1

)

!1!

)

y1...

((

m− r

)

!r!

)

yr. (17) Clearly, p

(

y0,...,yr |n,r,m

)

= h

(

y0,...,yr |n,r,m

)

 z0,...,zrh

(

z0,...,zr |n,r,m

)

,

wherethesuminthedenominatorextendsover allz0,...,zr

sat-isfying(12)–(14).Intheratioallfactorsonlydependingonn,r,m

cancelandwearriveat

Theorem1. p

(

y0,...,yr

|

n,r,m

)

= r  i=0



m i



yi /yi! _r i=0



m i



zi /zi! , (18)

where thesumin thedenominator runs overall valuesof z0,...,zr

satisfying(12)–(14).

Withhindsight,formula(18)isverynatural;theithterminthe numeratorreflects thenumber ofways ofordering i red ballsin eachofy_iurns,whereallurnscontainmballs.

LetAn,r,mbe thesetofall

(

r+1

)

-tuples

(

y0,...,yr

)

satisfying

(12)–(14). We denote the distribution on An,r,m withprobability

function(18)by

μ

(n,r,m).

3.2.DistributionoftheMarkovchain

Nowwecandetermine

P

(

Yj_=yj _{for j=0,...,p}

₎

_,

wherefor j=0,...,p the vector yj₌

₍

_yj

0,y

j

1,...,y

j

m_j

)

isan

arbi-trarypossiblevalueofYj₌

₍

_Yj

0,Y

j

1,...,Y

j m_j

)

.

TheMarkovchainclearlystartswithinitialdistribution

P

(

Y0_=l

₎

₌



n l



ε

l

₍

₁₋

_ε

₎

n−l_,

sincebeforewebeginpooling,thereisaninitialsupplyofnitems whereeach ofthemhas,independently oftheothers,probability

ε

ofbeingcontaminated.Theone-steptransitionprobabilities,i.e., theconditionaldistributions_PYj+1_|Yj,aregivenby

Theorem2. Letj_∈

{

0_,...,p_{− 1}

}

.Ifyj_andyj+1_arepossible

realiza-tionsofYj_andYj+1_{,respectively,then}

P

(

Yj+1=yj+1

|

Yj=yj

)

=



˜m_rj =1

μ

(

mj,r,mj+1

)

∗yj r



(

{

yj+1

}

)

. (19)

Here

μ

(n,r,m)∗k denotes thek-fold convolutionof

μ

(n,r,m) withitself,andtheright-handsideof(19)istheprobabilityofthe one-point set

{

yj+1

_}

_{under the convolution of the mj probability}

measures

μ(

mj,r,mj+1

)

∗y j r, r=1,...,mj. Here ˜ m_j r=1 denotes the convolutionproduct.

Proof of Theorem 2.. To have y_lj+1 groups in the

(

j+1

)

th cy-clethat haveexactlyl contaminateditems,y_lj+1 must bethe sum

ofthe numbers of subgroups with l contaminated items formed by splittingthe groupsafter the jth testing. At that time the re-maining contaminated groups are of size mj and are split into

subgroupsof size mj+1. There are yrj groups containingexactly r

contaminateditems, r=0,. . .,mj.The numbers ofnewly formed

subgroupswithexactlylcontaminateditemsfromallthesegroups haveto beadded toobtainthe numberofsuchsubgroupsin the

(

j+1

)

thcycle.Thesenumbersareindependentrandomvariables, whichexplainstheconvolutions.Eq.(19)isproved. _

The joint distribution of the finite nonhomogeneous Markov chaincannowbegivenexplicitly.

Theorem3. P

(

Yj_=yj _{for j=0,...,p}

₎

=



n y0



ε

y0

(

₁−

ε

)

n−y0



p_−1 j=0 ˜ mj r=1

μ

(

mj,r,mj+1

)

∗yj r



(

{

yj+1

}

)

. (20) 4. ExactdistributionsofthestoppingtimeandofM,N,N∗

The testing procedures are stopped when the subgroups are gettingtoosmall,sayafter

τ

1 cycles(this meansthat mτ1≥ ˆm>

m_τ1+1,cf.(3)).Weconsiderstoppingtimes

τ

=min

(τ

1,

τ

C

)

,as

de-fined in (1 )–(5). Now observe that the total number of tested groupsinthefirstjcyclesisgivenby

j  i=1 mi  l=0 Y_li. Consequently,cf.(4),

τ

C=sup

{

j

|

j  i=1 mi  l=0 Y_li≤ C

}

. (21)

Thusthedistributionof

τ

Ccanbeexpressedintermsofthe

under-lyingMarkovchain,whosedistributionwasobtainedinTheorem3. ThesameholdsforthedistributionsofM,NandN∗.Theexact for-mulasaregiveninthefollowingtheorem.

Theorem4. P

(

τ

C>j

)

=P



_j  i=1 mi  l=0 Yi l ≤ C



, (22) P

(

M=q

)

=τ 1−1  j=1 P



_j  i=1 mi  l=0 Yi l =q, j+1  i=1 mi  l=0 Yi l >C



+P



τ1  i=1 mi  l=0 Yi l =q



, q≤ C, (23) P

(

N=h

)

=τ 1−1  j=1 P



_j  i=1 miY0i=h, j  i=1 mi  l=0 Yi l ≤ C< j+1  i=1 mi  l=0 Yi l



+P



_τ 1  i=1 miY0i=h, τ1  i=1 mi  l=0 Yi l ≤ C



, (24) P

(

N∗=s

)

= n  z=s P

(

N=z

)



z s



(

1−

γ

)

s

_γ

z−s_{. (25)}

Proof.Eq.(22)follows immediatelyfrom(21).Next,M (thetotal numberofconductedgrouptests)canbedecomposedas

M=τ i=1 mi  l=0 Y_li, sincemi

l=0Yliisthetotalnumberoftestsconductedinthejth

cy-cle.Thisyields(23).Similarly,thenumberofitemstestedcleanat theELISAstationduringtheiterationscanberepresentedas N=τ

i=1

miY0i,

sothat itsdistributioncanalsobe obtainedfrom(20).Eq.(24)is noweasilychecked.Finally,(25)followsfromthelawoftotal prob-ability.

All theexpected valuesappearing inthe objectivefunction P˜,

whichwasintroducedin(9),cannowbeexpressedintermsofthe underlyingMarkovchain.Theresultingformulasareverylengthy. Theyarecomposedofcomplicatedconvolutionpowersof probabil-itymeasuresgiveninturnbyratios ofproducts ofbinomial coef-ficients,andareapparentlynotsuitablefornumericalcalculations onpresent-dayordinarycomputers. Thismayofcoursechangein thenearfuture,butinthemeantimethenaturalapproachisto re-sorttoapproximationswhichmaketheformulasmoreaccessible. Thisiscarriedoutintherestofthepaper. 

5. Anapproximationfornumericalpurposes

5.1. Somerealisticassumptionsandrelatedapproximations

Theexplicitformulasforthecomponentsoftheobjective func-tion,asderivedintheprevious section,are verycomplicated.For the purposes of optimization it is hence important to come up with approximations for these components. A starting point for such approximationsisthefollowing consideration.Thefollowing arerealisticnumbers(cf.Section1.3):

n=6720per week,



=P

(

unitcontaminated

)

=6.8× 10−4_,

γ

=4× 10−5_,

m1=48.

Sowehaveabout140initialbatchesperweek,and7280peryear. Theprobabilitythata48-batchisnotcontaminatedis

(

1−



)

m_1≈

1− m1



≈ 0.97.The probability that a 48-batch has one

contam-inated unit is m1

(

1−

)

m1−1≈ 0.03. Thus on average this

oc-curs approximately (only) four times per week. The probability thata48-batchhastwocontaminatedunitsis 1

2m1

(

m1− 1

)



2

(

1−



)

m1−2≈ 5× 10−4_{.Therefore,onaveragethisapproximatelyoccurs} onlyfourtimesperyear.Fromthisonegetsafeelingforthe pro-portions.

Now looking atblood units,instead ofbatches, per week, the numberofcontaminatedunitsinaweek isbinomiallydistributed with parameters n and



, and for n=6720 and



=6.8× 10−4

thisisextremely accurately representedby a Poissondistribution withparameter

λ

=n



=4.57(cf.Section6.5ofFeller,1968).This bringsustoour

Approximation Assumption1. Thedistributionofthenumberof contaminated itemsper week canbe approximatedby a Poisson distributionwithparameter

λ

,whichcanbeconsideredasthe ar-rivalrateofcontaminateditems(perweek),andwhichforgivenn

equals

λ

=n



.

So with X the number of contaminated items in an arbitrary week:

P

(

X=i

)

=e−λ

λ

i_{/i!,i=0,1,..., (26)}

Fig. 2. Poisson approximation versus Bernoulli simulation.

For ourobjective function (9) we need approximative expres-sionsfortheterms

E[N∗I₍N∗_≤d)], P

(

N∗>d

)

, E



(

N∗− d

)

I_(N∗_>d)



, E[M], E[N], E[

(

d− N∗

₎

_I(

N∗_≤d)]. (27)

Based on thefact (see above)that forthe realistic choice



= 6_.8_{× 10}−4 theprobability that a 48-batch containstwo contami-nateditemsapproximatelyequals5× 10−4_{,andthatinitialbatches}

arenotlargerthan64,weintroduce

Approximation Assumption 2. Onecan ignore the eventthat at leasttwo ofthecontaminateditemsareinthesameinitial batch (whichisaneventwithprobabilityoforder



2_).

We claim that both approximation assumptions are extremely accurate for realistic values of the contamination rate



. This is confirmedbyasimulationexperiment,requestedbyoneofthe ref-erees,inwhichwegenerated10,000weeks,eachwith140batches of size 48. In Fig. 2 we display the resulting distribution of the number of contaminated items per week, comparing it withthe Poisson approximation. In the simulation, 711 out of 1.4 million batchescontainedatleasttwocontaminateditems,supporting Ap-proximationAssumption2(thebinomialprobabilitiesinthe begin-ningofthissectiongiveaprobabilityof0.00051ofhavingatleast two contaminated itemsin abatch, leadingto an averageof 715 batchesper1.4millionwithatleasttwocontaminateditems).

So now assume that initially there are X=i contaminated items, allbelongingtodifferentbatches.Inthiscasetherewillbe

i contaminated batches inevery recycling (eachcontaminated by one item)so that, when there are h cycles,the number of tests in thisweek equals N_h= n

m₁+i

(

m1 m₂ +· · · + m_h−1 m_h

)

(remember that kj= m_j−1

m_j isthe numberof batchesbeing testedinthe jthcycle,

j₌2_,3_,...).Itfollowsfrom(2)and(4)thatthestoppingtimes

τ

1

and

τ

C,andthusalso

τ

,areconstants(dependingoniandonthe

batchsizes):wehave

τ

≡ hi,wherehi=min

(τ

1,

τ

C

)

.

LetusnowturntothedistributionofN,forgivenn.Firstofall,

P

(

N=n

)

=

(

1−



)

n_{≈ e}−λ_{; indeed,this is the casethat all items}

of this week are clean. Second, P

(

N=n− mh₁

)

equals the

prob-abilitythatthisweekthereisexactlyonecontaminateditem,i.e., P

(

X=1

)

.HenceP

(

N=n− mh₁

)

≈

λ

e−λ.Indeed,ifwetestthe

con-taminatedbatch

τ

times(insmallerandsmallerbatches),thenthe numberofitems(outoftheoriginal m1) that we donotsend to

PCRequalsm_τ.Generally,

P

(

N=n− imhi

)

≈ e

−λ

_λ

i_{/i!, i=0,1,... . (28)}

This corresponds to the probability of having X=i≥ 1 contami-nateditemsin thisweek; recallthat we ignoretheeventthat at leasttwoofthemareinthesameinitialbatch(havingprobability oforder



2_).

Let us next consider the distribution of M, the number of groupstestedinaweek.Withthesamereasoningasabove,again ignoringtheeventthatatleasttwocontaminateditemsappearin thesameinitialbatch,wehave:

P



M= n

m1 +

i

(

k2+· · · +khi

)



=e−λ

λ

i_{/i!, i=0,1,... . (29)}

Finally, the distribution of N∗ can be approximately determined from(25).Inview ofthefact that

γ

≈ 4× 10−5_{,wecan very}

ac-curatelyapproximate

P

(

N∗=s

)

≈ P

(

N=s

)

+

(

s+1

)

γ

P

(

N=s+1

)

, (30)

oreven_P

(

N∗₌s

)

_{≈ P}

(

N₌s

)

;infact,weproposetousethelatter formulain thefour terms inthe objectivefunction involving N∗. From (28)to (30)approximations forall terms in(27) are easily obtained. Thisyields a simple approach tothe objectivefunction thatcanbeusedforatentativeoptimization.

Thesamereasoningcanbeusedinthecasewhenthestopping timeisaprespecifiedconstant,sayh0_.Then

P

(

N=n− imh0

)

≈ e−λ

λ

i/i!, i=0,1,... . (31)

ThisagaincorrespondstotheprobabilityofhavingX=i contami-nateditemsinthisweek.Takingmeans,weget

E[N]≈ n−

λ

mh0. (32)

Furthermore,inthiscasewegetforthedistributionofM

P



M=_mn

1 +

i

(

k2+· · · +kh0

)



Inparticular, E[M]≈_mn 1 +

λ

(

k2+· · · +kh0

)

= n m1+

λ



_m 1 m2 +· · · + mh0₋₁ mh0



. (34)

5.2.Threepossibleapproximationcases

Consequently, taking into account Approximation

Assumptions 1 and 2 above, we may distinguish (at least)

threecases.

1.CaseI:thestoppingtimeisaprespecifiedconstanth0

Thisisanaturalandimportantcase.Weelaborateonthiscase inthesequel.

2. CaseII:X=iisgiven

Supposewedoafirstcyclewithgroupsofsizem1(adecision

variable). After this cycle, we count the number of contami-nated groups. Suppose this number, X, equals i. Ignoring the possibilityof havingmorethan one contaminateditem inthe samegroup orsubgroup,ifwecontinueuntilwe havedonea totalnumberofhicycles,weshallhaveN=n− imh_i.Our

deci-sionvariablesare m1,m2,...,mh_i (aswell ashi,unlesswe

de-cidethatwedoasmanycyclesaspossible,aslongasmh≥ ˆm

andNh≤ C).

Now, we claim that there are not so many cases to choose among.Ingeneral, thenumberofpossibilitiesdependsonthe prime factorization of n. The numbers of groups, i.e., the ki,

canbegeneratedbyanypartitionoftheprimefactorsofn, re-sultinginahugenumberofgroupingschemes.However,their number is drastically reduced by the upper bound on m1. In

our numerical examples we take n=6720 and suppose that

m1 ∈ {16, 20, 21, 24, 28, 32, 40, 42, 48, 56, 64}, that m2∈

{

8,10,12,...,m1/2

}

andthatm3∈

{

4,5,...,m2/2

}

.Wemay

al-low afourthandevenafifth cycle,butthat seemsnot realis-ticifm1 ≤ 64andwe makenewgroupsalways atleasttwice

assmallasthe groupsoftheprevious cycle. Furthermore,we cannot consider too small m1-values, since Nh ≤ C says that

n/m1 ≤ C.The smallnumberofcasesto be consideredmakes

it very easy to search among all possible cases. Searchingon theonehandmeans:checkingwhethertheconstraints regard-ingthestoppingtimesarenotviolated.Searchingontheother hand means: calculate the profit objective function, and take thelargestoneamongthoseforwhichthestoppingtime con-straintsarenotviolated.InCaseII,onehastoadaptthatprofit objectivefunction,givenbyformula(9),inanobviousway.We denoteitbyP˜itoemphasizeitsdependenceoni,andcalculate,

e.g.,itsfirsttermasfollows:

RwE[NI(N≤d)

|

X=i]=Rw

(

n− imhi

)

I(n−imh_i≤d).

In thisway we get a value forP˜i, forour giveni, andforall

combinationsofm1,m2,...,mh_iandhi.Asexplainedabove,one

cannowtakethecombination,amongtheadmissibleones,that yieldsthehighestprofit.

3. CaseIII:thegeneralcase

One way to treat this case is to take case II (but with hi

replaced by a decision variable H that does not depend on the actual value of X₌i), and multiply P˜_i by P

(

X₌i

)

and sumover i=0,1,...,z. Here zmight be such thatP

(

X>z

)

= ∞

i=z+1e−λ

λ

i/i!≤ 10−4,with,e.g.,

λ

=4.57.Dothisagainforall

possiblecombinationsofm1,m2,...,mH andnumberofcycles

H,whereweallowHtotakethevalues1,2,3,4,5,say,andwhere wecheckwhichcombinationsdonotviolateourstopping con-ditions.

Below we furtherelaborate onCases IandIII(Case II maybe viewedasaspecialcaseofCaseIII).

ElaborationonCaseI

Forsimplicity,we assume that n isdivisible by m_h0.All func-tionalsappearingin(11) canbesimplyobtainedbyusingthe ap-proximationsin(32)and(34).Indeed,byletting

˜ k=. n− d mh0 andF

(

t

)

=P

(

X≤ t

)

= [t]  i=0 e−λ

λ

i_/i!,

where[t]isthelargestnonnegativeintegerlessthanorequaltot

(recallthatXhasaPoissondistributionwithmean

λ

),the expres-sionsforP1,...,P6areobtainedasfollows.First,

P1=E[RwN∗I(N∗_≤d)] E[RwNI(N≤d)]=RwE[NI(N≤d)],

andfortherighthandsidenotethatE[NI₍N≤d)]=E[N]− E[NIN>d)],

where E[N]=n−

λ

m_h0 and E[NIN>d)]=nP

(

X=0

)

+

(

n−

m_h0

)

P

(

X=1

)

+...

Astraightforwardcomputationyields E[NI ₍N≤d)] =

⎧

⎪

⎪

⎨

⎪

⎪

⎩

n₋

λ

m_h0, ifd≥ n(k˜≤ 0) n(1− F(0))−

λ

mh0, ifd <nand0<k˜≤ 1 n



1− F(˜_{k− 1)}



₋

_λ

_m h0



1− F(_k˜_{− 2)}



_{, ifd <nand2≤ ˜k∈N} n



1− F([˜k])



−

λ

mh0



1− F([k˜]− 1)



, ifd <nand1<_k˜_∈/N. Next, P2=E



{

dRw+Rˆw

(

N∗− d

)

}

I(N∗_>d)



dRwP

(

N>d

)

− ˆRwE[

(

d− N

)

I(N>d)], where dRwP

(

N>d

)

⎧

⎨

⎩

0, ifd≥ n

(

_k˜_{≤ 0}

₎

dRwF

(

k˜− 1

)

, ifd <nandk˜∈N dRwF

(

[k˜]

)

, ifd <nandk˜∈/N, and E[

(

d− N

)

I₍N>d)]=E

(

d− N

)

− E[

(

d− N

)

I(N≤d)] =d−

(

n−

λ

mh0

)

− E[

(

d− N

)

I_(N≤d)].

Combiningtheaboveresultswefindthat E[(d− N)I ₍N≤d)]

⎧

⎪

⎪

⎨

⎪

⎪

⎩

d_−n+

λ

m_h0, ifd≥ n(k˜≤ 0) (d−n)(1−F(0))+

λ

mh0, ifd<nand0<k˜≤ 1 (d_−n)



1_−F(k˜−1)



₊

λ

m_h0



1_−F(k˜−2)



_, ifd_<nand2_{≤ ˜}k_∈N (d−n)



1−F([k˜])



+

λ

mh0



1−F([k˜]−1)



, ifd<nand1<_k˜_∈/N and E[

(

d− N

)

I₍N>d)

⎧

⎪

⎨

⎪

⎩

0, ifd≥ n

(

_k˜_{≤ 0}

₎

(

d− n

)

F

(

0

)

, ifd <nand0<_k˜_{≤ 1}

(

d− n

)

F

(

k˜− 1

)

+

λ

mh0F

(

k˜− 2

)

, ifd <nand2≤ ˜k∈N

(

d− n

)

F

(

[k˜]

)

+

λ

mh0F

(

[k˜]− 1

)

, ifd <nand1<k˜∈/N implyingthat P2 ⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩ 0, ifd≥ n(_k˜_{≤ 0)} dRwF(0)+Rˆw(n− d)F(0), ifd<nand0<k˜≤ 1 dRwF(k˜−1)+Rˆw[(n−d)F(k˜−1)−λmh0F(k˜−2)], ifd<nand2≤ ˜k∈N dRwF([˜k])+Rˆw[(n−d)F([k˜])−λmh0F([˜k]−1)], ifd<nand1<k˜∈/N. Theothercomponentsoftheprofitfunctionaregivenby

P3=Cpurn P4=CelisaE[M]=Celisa



_n m1+

λ



_m 1 m2 +...+ mh0₋₁ mh0



P5=CpcrE[N]=Cpcr

(

n−

λ

mh0

)

P6=Cpenalty

(

E[

(

d− N∗

)

I(N∗_≤d)]

)

Cpenalty

(

E[

(

d− N

)

I(N≤d)]

)

.

ElaborationonCaseIII

Forconvenience,anypossibleassignmentforthevaluesofthe– nownotprespecified– stoppingtimeh0_,m

1,...,mh0willbecalled a strategy,denotedby S₌S

(

h0_,m

1,...,mh0

)

.Foranygiveninitial valuesofparametersinvolved(n,Rw,Rˆw,...),thenumberof

strate-giesisfinite.ThesetofpossiblestrategiesisdenotedbyS.

Foreverystrategy Stheexpectedvalue oftherelatedprofitis givenby ˜ P

(

S

)

= z  i=0 P

(

X=i

)

× Pi,S, (35)

where X, as defined above, is the random variable counting the number ofcontaminated blood units(X _{∼ Poiss}(

λ

)) andthe con-stantzisdeterminedbytheconstraint

P

(

X>z

)

= ∞

i=z+1

e−λ

λ

i_{/i!≤ 10}−4_{, (36)}

and P˜i,S is the conditional expected profitwhen using strategy S

giventhatX=i.

Beforecontinuingourderivations wemakethefollowing com-ments:

1. Foranygivenvaluesoftheparameters involved,there arenot so many strategies. For example, for n₌6720_, m1 ≤ 64 and

minimalbatchsizemˆ=4,therearenomorethan99strategies. Infact,wewillprovideallofthemattheendofSection6. 2. As claimed in Section 5.1, we shall ignore the possibility of

having more than one contaminated item in the same group orsubgroup.IfX=i andwe continuerecyclinguntil we have done a total numberof hi cycles, we then shallhave N=n−

im_h i and M=Mi=



n m1 + i



m1 m2+· · · + mhi−1 mhi



. (37)

Ourdecisionvariablesarem1,m2,...,mhi (aswell ashi,unless

wedecidethatwedoasmanycyclesaspossible,aslongasmh≥

ˆ

m=4andNh ≤ C).

The smallnumber ofcases tobe considered makes it easy to searchamongallpossiblecases.Searchingontheonehandmeans: checkingwhethertheconstraintsregardingthestoppingtimesare not violated. Searching on the other hand means: calculate the profitobjectivefunction,andtakethelargestoneamongthosefor whichthestoppingtimeconstraintsarenotviolated.

Accordingly, one has to adapt the profit objective function in

(9)inanobviousway.Indeed,foragivenX=i,theexpectedprofit ˜ Piis ˜ Pi=Pi,1+Pi,2− Pi,3− Pi,4− Pi,5− Pi,6, where P_i,1 Rw

(

n− imhi

)

I(n−imh_i≤d)

(asX=iisgiven,theexpectationturnsouttobeaconstant),or P_i,1



Rw

(

n− imhi

)

, ifn− imhi≤ d 0, otherwise , Pi,2

(

Rwd+Rˆw

(

n− imhi− d

))

I(n−imhi>d), or Pi,2



Rwd+Rˆw

(

n− imhi− d

)

, ifn− imhi>d 0, otherwise , Pi,3=Cpur× n, Pi,4=Celisa× Mi, whereMiisgivenby(37), Pi,5=Cpcr

(

n− imhi

)

, Pi,6 Cpenalty



d−

(

n− imhi

))

I(n−imhi≤d)



or P_i,6



Cpenalty

(

d−

(

n− imhi

))

, ifn− imhi≤ d 0, otherwise .

6. Numericalandsensitivityanalysis

Inthissectionwe presentanumericalandsensitivityanalysis forCases IandIII. Since manydifferentparameters are involved, onecouldconsiderawiderangeofparametervalues,studyingthe influenceofeachofthemontheprofitfunction.Wemainlyrestrict ourselvestoparticularparametervalueswhichseemtoberealistic inthecaseoftheIsraeliCentralBloodBank.Thisallowsustofocus onafewkeyaspects,viz.:

(i)ForCaseI,we studytheinfluenceofdemand onprofit.For this case, we also consider the six components P1,...,P6 of the

profitfunctiontogetanimpressionoftheirrelativecontributions. (ii)InCaseIIIwestudytheeffectofdoingmultipletestsatthe ELISAstationontheprofitfunction.CaseIIIalsoallowsustostudy theeffectofgroupsizesattheELISAstationontheprofitfunction.

6.1. CaseI

Wefirstplotprofitversusdemandforh0_{=4andforfourcases}

ofiteratedbatchsizes,givenbythefollowingtable:

Case m1 m2 m3 m4

A 48 24 12 6

B 60 30 15 5

C 64 32 16 8

D 64 16 8 4

Demandchangesfrom5000to7000andthechosenparameter valuesare:

n=6720;

ε

=0.00068;Cpur=180; Celisa=2.5;

Cpcr=40; Cpenalty=25.

WeconsidertwocasesforRw andRˆw withRw=Rˆw=250and

Rw=250>Rˆw=230.

Byplottingprofitvs. demand forthe twocases ofRw andRˆw

one realizes that the behavior of the four different strategies is quitesimilar andtheprofit functionsdonot intersect.Indeed,as canbeseenfromFigs.3and 4,theprofitfunctionsbehave simi-larlyforeachofthefourselectedstrategies,andstrategyDhereis thebestforalldemandvalues.IfRw=Rˆw,thentheprofitfunction

remainsalmostconstantuntiln≈ d;forlargervaluesofd,i.e., de-mandcannot be satisfied, the penalty factorP6 causesthe profit

function todecrease linearly.IfRw>Rˆw, then theprofit function

increaseslinearlyinduntiln≈ d;thisisduetotermP2.Forlarger

valuesofd,factorP6 againcauses theprofitfunctionto decrease

linearly.

Variousothervaluesoftheparameters thathavebeen consid-ered show a pattern similar to the one presented. In particular, whenRw=Rˆw andCpenalty=0,theprofitstaysconstantasa

func-tionofthedemand(indeed,noticethatnowP1+P2=RwE[N]does

notdependond,andneitherdoP3,...,P6).

Next,westudytheimpactonthecostandpenaltycomponents oftheexpectedprofitP˜whenslightlyloweringthedemand.

Theparametervaluesconsideredare:

ε

=0.00068;Cpur=180;Celisa=2.5;Cpcr=40;