Screening and merging algorithm for the detection of copy-number alterations

0.14 MB·

by  Murilo S. Pinheiro

BiometricalJournal52(2010)61,zzz–zzz/DOI:10.1002/bimj.200100000
Screening and merging algorithm for the detection of copy-
number alterations.
MuriloS.Pinheiro,BeniltonS.Carvalho,andAlu´ısioS.Pinheiro∗
7 DepartofStatistics,IMECCUniversityofCampinas,Se´rgioBuarquedeHolanda,65113083-855. Camp-
1 inasSPBrazil.
0
2
Receivedzzz,revisedzzz,acceptedzzz

n
a
J Wecallchange-pointproblem(CPP)theidentificationofchangesintheprobabilisticbehaviorofasequence
ofobservations. SolvingtheCPPinvolvesdetectingthenumberandpositionofsuchchanges. Ingenetics
6
thestudyofhowandwhatcharacteristicsofaindividual’sgeneticcontentmightcontributetotheoccurrence
1
andevolutionofcancerhasfundamentalimportanceinthediagnosisandtreatmentofsuchdiseasesandcan

be formulated in the framework of chage-point analysis. In this article we propose a modification to a
]
P existingmethodofsegmentationwiththeobjectiveofproducingaalgorithmthatisrobusttoavarietyof
A samplingdistributionsandthatisadequateformorerecentmethodofaccessingDNAcopy-numberwhich
mightrequirearestrictionontheminimumlengthofaalteredsegment.
.
t
a
t Keywords: SNPArray;Change-Point;Copy-number;SaRa;segmentation;
s
[

1
v
1 Introduction to Change-PointProblems and their RoleinGenetics
5
0
4 Theproblemofidentifyingachangeintheprobabilisticbehaviorofasequenceofobservationsisusually
4 called the change-pointproblem (CPP). Solving the CPP involves detecting the number and position of
0 suchchanges. Itiscommontoassumethatthechange-pointstructureissuchthatstochasticpropertiesof
.
1 thesequenceareapproximatelypiecewiseconstant.
0 Let{X ,...,X }beasequenceofindependentrandomvariables(i.r.v.).Wesaythat{X ,...,X }is
1 n 1 n
7
segmentedwithrespecttothesequenceofchange-points1< τ < ...<τ ≤nif,andonlyif,there
1 1,n s,n
ares+1cumulativedistributionfunctionsF ,F ,...,F suchthatF 6= F , forj ∈ {0,...,s−1},
: 0 1 s j j+1
v andthecumulativedensityfunctionofeachX ,F ,canbewrittenas:
i Xi
i
X
F if1≤i<τ
r 0 1,n
a  F1 ifτ1,n ≤i<τ2,n
FXi = ...
F ifτ ≤i≤n
 s s,n
Inthiscontextwesaythatτ ,...,τ arethechange-pointsin{X ,...,X }.
1 s 1 n
Thepreviousgeneralformulationadmitsthreeimportantrestrictions: 1)iftheF ,...,F areassumed
1 s
to belongto the same location family then the associated CPP is called the mean change-pointproblem
(mCPP);2)ifF ,...,F belongtothesamescalefamilythenwehavethecalledvariancechange-point
1 s
problem(vCPP)and;3)ifthecumulativedistributionfunctionsbelongtothesamelocation-scalefamily
thenwehavethecalledmeanandvariancechange-pointproblem(mvCPP).HerewefocusonthemCPP.
Assuming that s is the number of change-points in a sequence the statistical problem of testing the
existenceofatleastonechangeisequivalenttothefollowinghypothesistest: H : s = 0;H : s > 0. If
0 1
∗Correspondingauthor:e-mail:[email protected],Phone:[phone],Fax:[phone]
(cid:13)c 2010WILEY-VCHVerlagGmbH&Co.KGaA,Weinheim www.biometrical-journal.com
2 Pinheiroetal.anddd:Runningtitle
weareunabletorejectH thannothingelsehastobedoneexceptfindestimatesforthedesiredquantities;
0
otherwise,findinggoodestimatesforthechange-pointsisanecessarystepforanydownstreamanalysis.
Solutions to the CPP are of great importance in many areas of research such as climatology, finance
and genetics and, for that reason, a great variety of solution methods exist. The study of change-point
problems started in Page (1954) where a approach based on cumulative sums (CUSUM), further devel-
oped in Page (1955), was proposed. After this initial study the change-point problem received a large
amountofattentionresultinginnumerousmethodsandapproaches,suchas:useofthemaximaofstatistics
calculatedalongthesequenceofobservations(VenkatramanandOlshen, 2007;KirchandMuhsal, 2014;
Fryzlewiczetal., 2014; BrodskyandDarkhovsky, 1993; Fryzlewiczetal., 2014), optimization methods
(AugerandLawrence,1989;Killicketal.,2012)andhiddenMarkovmodels(AugerandLawrence,1989;
Killicketal., 2012; Yau, 2013; Wangetal., 2007; Colellaetal., 2007; Beroukhimetal., 2006; Haetal.,
2012).Thoseexamplesdonotfullyencompasstheextentofchange-pointanalysisresearch.KhodadadiandAsgharian
(2008)providesaextensivereviewonvarioussubfieldsofchange-pointanalysis.
Therearemanyareasofappliedresearchthatbenefitsfromchange-pointanalysisethnics.Instatistical
processescontroloneisinterestindetecting,asfastaspossible,thepointafterwhichaprocessesstopped
behavingintheexpectedway(Montgomery,2009). Itisclearthatchange-pointanalysistechnicscanbe
appliedindetectingthistypeofchange(Bassevilleetal.,1993). Otherareawherechange-pointanalysis
playsa role iseconometricswhere itis importantto detectchangesin the behaviorof time ordereddata
(BaiandPerron,2003;Allenetal.,2013).
Ingeneticsthestudyofhowandwhatcharacteristicsofaindividual’sgeneticcontentmightcontribute
to the occurrenceand evolutionof cancer has fundamentalimportancein the diagnosisand treatmentof
suchdiseases. Ofparticularinterestinthiscontextarethecopy-numberstructureofagenome.
Inmedicineresearch,speciallygenetics,change-pointanalysisalseplaysacrucialroleintheunderstand
of complexdiseases (diseases that are consequenceof a complexinterationsbetween multiple genetical
andenvironmentalfactors). Theinventionofarraytechnologiesenabledtheanalysisofalargenumberof
DNAlocus(witharrayCGHouarraySNP)oreventhewholegenome(withwholegenomesequencing).
Usually the results of a array experimentare two numericalsequences one of those related to the DNA
copy-numberstructure.Changesinthissequencelevelindicatesincrementsordecrementsintheanalyzed
DNA genomic content. Those changes have been shown to effect the occurrence and development of
cancers.
Theevolutionoftechnologicalcapabilitiesmakesthesegmentationofsuchsequencesevenmorechal-
lenging since modern methods generate data with a large numbers of observations, possibly more than
3billion observations(one for each base in the DNA chain). Searching for change-pointson such large
data sets can only be feasible with methodswhose computationalcomplexityare as close as possible to
linear. Anotherissue is thatasthe technologicalevolutionincreasesthedetectionof smallsegmentsbe-
comesless practicallyrelevant, Mason-Suaresetal. (2013), andmightproduceunreliableresults since a
biggerfalsepositiverateisinevitableifthechosenmethoddonotpreventthedetectionofsmallsegments.
We proposean adaptionto a existingmethod,(NiuandZhang, 2012), to the scenariowherethelarge
amountofdatahinderstheapplicationofotherslowermethods,wheretheadoptionofaminimumlength
for the distance between segments is reasonable and the sample distribution deviations from normality
hinderstheapplicationofmethodsrelyingonthenormaldistribution.
{X ,...,X } is a sequence of independent random variables and we take µ = E(X ) and σ2 =
1 n i i
V(X ) < ∞. We also adopt the following notation: [i : j] = {i,i+1,...,j} if i and j are natural
i
numberssuchthati<j. IfI ⊂{1,...,n}then
1
X = X ,
I i
|I|
Xi∈I
where|I|denotesthecardinalityofI.
(cid:13)c 2010WILEY-VCHVerlagGmbH&Co.KGaA,Weinheim www.biometrical-journal.com
BiometricalJournal52(2010)61 3
2 The SaRa Segmentation Procedures
NiuandZhang(2012)proposedamethodbasedontheintuitionthatthepropertyofaparticularposition
i, alongthesequence{X ,...,X }, beingachange-pointcanbe assessed byconsultinga vicinityof i.
1 n
Basedonthissimpleobservationtheauthorsusethecumulativesumofneighboringregions:
i−1 i+k−1
Lk =(cid:12) X − X (cid:12),
i (cid:12) j j(cid:12)
(cid:12)j=Xi−k Xj=i (cid:12)
(cid:12) (cid:12)
(cid:12) (cid:12)
wherei = k+1,...,n−k+1toquantify(cid:12)theevidenceoftheex(cid:12)istenceofachange-pointinpositioni.
WedefineV(i)=[i−k,i+k−1]fornotationeasy.
Theprocedureforcallingchange-pointsisbasedontwosteps:
• Screening step; change-pointcandidatesare selected by finding the valuesof i such that Lk ≥ Lk
i j
forallj ∈ V(i)andLk > δ whereδ isapredefinedrealconstant. Thatis,themethodselectallthe
i
globalmaximainthevicinitiesV(i)whicharegreaterthanathresholdδ. Thevalueofδissuggested
toberelatedtotheminimalvalueofashiftinthemean.WewilldenotebyJ thesetofallcandidate
change-pointsfoundinthisstep.
• Rankingstep;Thenextstepistherankingofallcandidatechange-pointsonJ withrespecttosome
information criteria. In NiuandZhang (2012) the authors adopt the Bayesian information criteria
(BIC)andthefollowinggreedheuristictoselectafinalsetofchange-pointsisproposed. LetJ be
(i)
the resultofremovingthe i-thcandidatechange-pointfromJ. If noset J hasa lowerBIC than
(i)
J thenJ isthefinalsetofchange-points. OtherwiseletibetheindexsuchthatJ isthesetwith
(i)
minimalBICvalue,thenJ isupdatedtoJ andwerestartthesecondstepwiththisreducedsetof
(i)
change-pointscandidates.
TherearesomeissuesinutilizingtheSaRaalgorithmasproposedbyNiuandZhang(2012). Amean-
ingfulvalueofδ mightnotbeeasilyselectedbecauseimpuritiesinthecellsamplesfromwhichtheLRR
sequencesweregeneratedchangetheminimaldifferencebetweenneighboringsegmentsmeans. Alsothe
choiceofBICasacriteriaforremovingincorrectelyidentifiedchange-pointcandidatesmakethemethod
model dependent and deviations from the normal distributions might result in unwanted results. Other
problemswiththeBICcriteriawereidentifiedinXiaoetal.(2014).
Xiaoetal. (2014) proposeimprovementsto overcomethese difficulties. The authorsproposethe use
of quantilenormalizationin a prepossessingstep to addressdeparturesfromnormality. Theselection of
δ isdoneviaaMonteCarlosimulationprocedureforfindingaparticularquantileofthemaximumofLk
i
distribution.Theauthorsalsonotedthattherankingproceduredidnotremovefalsepositivesinasatisfac-
toryway.Theythenproposedanewprocedureforremovingincorrectlydetectedchange-points.Thisnew
approachusestheclusteringinformationofallsegmentsgivenbytheinitialchange-pointscandidatesJ
in,atmost,threegroupswiththehelpofthemodifiedbayesiancriteria. Whentwoneighboringsegments
areassignedtothesameclusterthechange-pointbetweenthemisremovedfromJ.
ThemodificationsformulatedbyXiaoetal.(2014)resultedinastrategythatismorerobusttodeviations
from normality than the SaRa method. One disadvantage, noted by the authors, is that the number of
clusterschosentogroupthesegmentsmayaffectthesegmentationwhenthetruenumberofclusterdiffer
fromthechosenvalue. Anotherissue isthatthe use ofquantilenormalizationmayundulyinfluencethe
analysisofmicroarrayexperimentsdata,asnotedbyQiuetal.(2013). Anotherdifficultyistheexecution
timeburdenofperformingMonteCarlomethodsforlargesamples.
Inthenextsectionweproposeastrategytoselectthethresholdδwhenwedonothaveanyinformation
regarding the minimum mean change required to identify a change-point. We also propose a novel ap-
proach,robusttosampledistributionandsize,toreducethenumberofincorrectlydetectedchange-points.
Ourproceduresarebestsuitedforsituationswheretheminimallengthofasegmentisnotsosmall(atleast
(cid:13)c 2010WILEY-VCHVerlagGmbH&Co.KGaA,Weinheim www.biometrical-journal.com
4 Pinheiroetal.anddd:Runningtitle
25observations).Thelengthrequirementisasugestionmadebymicroarraymanufacturers,likeAffymet-
rics,toreducefalsepositivesondownstreamanalysis.Despitetherequirement,thisisacommonsituation,
giventhatrecentmicroarrayproductscontainmillionsofprobes(Mason-Suaresetal.,2013).
3 Screening and Merging Segmentation Procedure
WepresentamethodcloselyrelatedtotheSaRamethod.Ourmodificationsare:themeasureusedtoassess
whether or not a position is a change-points, the choice of threshold δ and the procedure for removing
incorrectlyselectchange-pointscandidates.
3.1 ThresholdSelectionandScreeningProcedure
WemodifyLk bydefiningthefollowingstatistic:
i
1 2 −1/2
Mk = X −X ,
i S (cid:18)k(cid:19) [i−k:i−1] [i:i+k−1]
seq (cid:12) (cid:12)
(cid:12) (cid:12)
where
n
1
S2 = (X −X )2,
seq 2(n−1) i i−1
Xi=2
isaconsistentestimatorofσ2,eveninthepresenceofmultiplechange-points,givenafewweakconditions
onthechange-pointsevolutionasnincreases.
WecarriedthescreeningstepofSaRainthesamewaybutreplacingthestatisticLkwithMk. Itallows
i i
asimple,intuitiveandfastwayofselectingathresholdδ.
IfweassumethattheX arenormallydistributed(orthatkisbigenough)andthattherearenochange-
i
points along the sequence, then each Mk is approximately distributed as a folded normal distribution,
i
havingcumulativedistributionfunctionF givenby
fold
F (x)= x 2e−u2/2d,u (1)
fold Z rπ
0
WeproposetheuseoftheuppertailquantilesofF ,i.e.,forα∈(0,1)takeδ(α)=F−1 (1−α)as
fold fold
thethresholdforthescreening.
AspointedbyNiuandZhang(2012),thescreeningprocedurecanbewithmorethanonevalueofkand
wealsoutilizethisapproachinourmethod.
3.2 MergingProcedure
Afteraseriesofscreeningproceduresisappliedto{X ,...,X },weareleftwithasetofchange-point
1 n
candidates,J ={τˆ ,...,τˆ },whichmightcontainchange-pointsthatdonotsatisfytheconditiononthe
1 p
minimumsegmentlengthandcandidatesthatrepresentfalse-positivedetections.Tofilterthosecasesfrom
J we adopt the following procedure: we sequentially test if the mean of observations on the left of τˆ
i
differsfromthemeanofobservationsontherightofτˆ. Thatis,foreachiwecomputetheusualstatistic
i
forcomparingthemeanoftwogroups:
−1/2
1 1 1
T = + (µˆ −µˆ )
S (cid:18)τˆ −τˆ +1 τˆ −τˆ (cid:19) [τˆi−1:τˆi−1] [τˆi:τˆi+1−1]
seq i i−1 i+1 i
andcompareittothecriticalvaluefoundbythenormalapproximationtothedistributionofT.
(cid:13)c 2010WILEY-VCHVerlagGmbH&Co.KGaA,Weinheim www.biometrical-journal.com
BiometricalJournal52(2010)61 5
If T fall outside a predefinedconfidenceintervalthen we proceedto the nextchange-pointcandidate
τˆ . Otherwiseweremoveτˆ fromthelistofchange-pointscandidatesandsincethischangetheobserva-
i+1 i
tionstotherightofτˆ werepositionτˆ with
i−1 i−1
−1/2
1 1
τˆ = argmax + µˆ −µˆ ,
i τi−1+k′≤j≤τi+1−k′+1(cid:18)j−τˆi−2 τˆi−j(cid:19) (cid:12) [τˆi−2:τˆi−1−1] [τˆi−1:τˆi−1](cid:12)
(cid:12) (cid:12)
whichistheestimatorofasinglechange-pointbasedonthelikelihoodratiostatisticforthecaseofnormally
distributeddata(SenandSrivastava,1975). Notethatthevalueofk′doesnotneedequalk.
Becausethismethodsisbasedonmergingadjacentsegmentsthatfailtorejectthehypothesisthatthey
share the same mean we call this step the merge step. And, for the same reason, we call our method
screeningandmerging(SaMe)segmentationmethod.
4 SimulationStudies
Wenowturntotheissueofassessingtheabilityofourproposedmethodstocorrectlyidentifythelocation
ofchange-pointsinthemeanalongasequenceofobservations.Alsoweareinterestedinhowthestochastic
distributionofsaidobservationsmightaffectthemethod’sperformance.
WecompareourmethodSaMetothewellestablishedproceduresproposedinVenkatramanandOlshen
(2007),whichwedenotebyCBS,andKillicketal.(2012),whichwedenotebyPELT.
Theparameterchosenforourmethodareα =0.01bothinthescreeningstepandinthemergingstep.
Weutilizeamultiplebandwidthsapproachwithk ∈{25,50,100}inthescreeningstep. Thevalueofk′is
themergestepissetto20.
In the application of CBS we followed the authorsrecommendationsadopting a level of significance
of0.01. TheparameterschosenforPELTwilldependontheparticulardatadistribution. Inallsituations
we adopta minimumdistances between change-pointsof 20 (for consistency with the value of k′). For
normally distributed data we chose the BIC penalization (since choosing AIC would result in too large
false-positivesdetection)andinthecasewherewesampleobservationsfromarealdatasetwechosethe
AICpenalization(sinceinthiscasetheBICalternativeresultinaverylowrateofchange-pointdetection).
Thecriteriaforaccessinghowfastamethodestimatethechange-pointswasitsexecutiontimeinsec-
onds(t). The accuracyof a methodis accessed in the followingway, if τ is the indexassociated with a
change-pointwe consider two windows containing τ: Wtol = {τ −tol,...,τ,...,τ +tol} where tol
τ
is ether10or 5. We say thata methoddetectedτ with tolerancetol if thereis a estimatedchange-point
insideWtol. Wedenotebyp10thepercentageofchange-pointsthatwheredetectedwithtolerance10by
τ
themethod.Alsop5istheanalogouswithtolerance5. ThelastcriteriaisFPdenotingthemeandifference
betweentheactualnumberofchange-pointsandthenumberofrealchange-pointsdetectedwithtolerance
10. Thislastcriteriaisawayofaccessingtherateoffalsediscoveriesofchange-points.
4.1 NormalyDistributedData
Inthefirstsetofsimulationsweconsiderasequenceofobservations{X ,...,X }withn=10,000and
1 n
suchthat
X =f(i)+ǫ ,
i i
where the ǫ are all N(0,1) and the function f specify the change-pointstructure. To fully investigate
i
theimpactofchange-pointspositionsontheresultsofthevariousmethodstestedtheformoff foreach
particular simulation is chosen in the following way. First the number of change-points m along the
sequence is fixed to be one of 2, 4 or 6, the length of the segmentdetermined by the change-pointsl is
chosentobeonofthevalues25,50or100,finallythemeanshiftsisfixedatonethevalues1,1.5or2.
Thenwerandomlysamplem/2naturalnumbers,cp ,...,cp ,uniformlyform[(i−i)L+l:i∗L−l),
1 m/2
where L is the greatest integer smaller than 10000/(m/2). The change-points positions are chosen as
(cid:13)c 2010WILEY-VCHVerlagGmbH&Co.KGaA,Weinheim www.biometrical-journal.com
6 Pinheiroetal.anddd:Runningtitle
{cp ,cp +l,...,cp ,cp +l} and f(i) is given by f(x) = s, if cp ≤ i < cp +l, for some
1 1 m/2 m/2 i i
i ∈ {1,...,m/2}, and f(x) = 0 elsewhere. For all the 27 possible combinations of m, l and s we
construct100randomsamples.
Thefirstthingtonoteisthatallmethodsperformwellintermsofp5andp10whens = 1.5ors = 2,
and for this reasons those results are omitted. It is not surprising since the presence of a change-point
becomesmoreevidentthegreatersis. Forthevaluesofs=1therearedifferencesbetweenthemethods.
Table1 ComparisonofSaMe,CBSandPELTonnormallydistributeddatawithchangesinmeanofmag-
nitude1. Intermsofexecutiontimeandchange-pointdetectionSaMeoutperformedbothCBSandPELT.
SaMeandCBShavesimilarfalse-discoveryratesinthiscase.PELThavethesmallestfalse-discoveryrate.
SaMe CBS PELT
t 0.1532 0.8936 0.4226
p10 0.5400 0.4500 0.2700
single
p5 0.4900 0.3700 0.2400
FP 0.0700 0.0600 0.0200
t 0.1386 1.0853 0.3106
p10 0.5500 0.5125 0.2675
double
p5 0.5050 0.4625 0.2575
FP 0.0900 0.1100 0.0500
t 0.1317 1.4698 0.2352
p10 0.5100 0.5183 0.2317
trip
p5 0.4717 0.4633 0.2167
FP 0.1100 0.2900 0.0500
Table 2 Comparison of SaMe, CBS and PELT on normally distributed data with changes in mean of
magnitude1. SaMehasthesmallestexecutiontime. Intermsofchange-pointdetectionCBSoutperforms
SaMeandPELTinallinstances. Falsediscoveriesratesofallmethodsaresimilar.
SaMe CBS PELT
t 0.1400 0.6906 0.3823
p10 0.8100 0.9200 0.8200
single
p5 0.7250 0.7950 0.7050
FP 0.2500 0.1600 0.1200
t 0.1330 0.8445 0.2669
p10 0.8900 0.9150 0.8500
double
p5 0.7975 0.8050 0.7450
FP 0.2000 0.4100 0.2600
t 0.1105 0.7965 0.2042
p10 0.8717 0.9367 0.8717
trip
p5 0.7850 0.8417 0.7600
FP 0.3300 0.4200 0.3200
(cid:13)c 2010WILEY-VCHVerlagGmbH&Co.KGaA,Weinheim www.biometrical-journal.com
BiometricalJournal52(2010)61 7
Table 3 Comparison of SaMe, CBS and PELT on normally distributed data with changes in mean of
magnitude1. Allmethodshavesimilarchange-pointdetectionabilityandfalse-positiverates. SaMeonce
againhasthesmallestexecutiontime.
size=100
t 0.1372 0.5050 0.3819
p10 0.8750 0.9250 0.9500
single
p5 0.7450 0.8450 0.8300
FP 0.3800 0.1900 0.1000
t 0.1195 0.5420 0.2668
p10 0.9300 0.9350 0.9400
double
p5 0.7875 0.8225 0.8450
FP 0.3300 0.3500 0.2400
t 0.1019 0.5903 0.1918
p10 0.9383 0.9600 0.9367
trip
p5 0.8200 0.8500 0.8383
FP 0.4500 0.4400 0.3900
Table1showsthatforthecaseofasinglealteredsegmentsofsize25andmeanshiftofmagnitude1the
SaMemethodoutperformallothermethodsintermsofthenumberofcorrectelyidentifiedchange-points.
When there are two or three altered segments SaMe still is the best methods in terms of change-point
identification ability, but CBS performs very close to it. PELT had a poor change-point identification
performancefor everynumberof altered segments. In terms of executiontime CBS took the longestto
completeexecution. SaMewasbyfarthefastestmethod. FalsediscoveryratesweresimilarforCBSand
SaMewhilePELTproducedthelestnumberoffalse-discoveries.
InTable2weseetheresultsforthecaseofalteredsegmentsoflength50andmeanshiftdifferenceof
magnitude 1. In this case we can see a slight increase in false-discoveryrates for CBS and SaMe. The
executiontimefollowthesamepatterwithSaMebeenthefastest. IntermsofdetectionabilityCBSwas
themostaccuratemethod.SaMeandPELThadasimilardetectionability.
Table 3 show the results for the case of altered segments of length 100 and mean shift difference of
magnitude1. TheresultsherearesimilartotheonesobtainedinTable2butnowthedetectionabilitiesof
eachmethodsareclosertooneanother.
4.2 SamplingRealData
We perform a second simulations study since the critical values for SaMe are approximationsbased on
the normaldistribution. Thisallow usto evaluatethe effectsofdeviationfromnormality. We performa
simulation study using observationssampled from a pool of data arising from a experimentin genetics.
ThedataisavailableintheacnrRpackage(Pierre-JeanandNeuvial,2016).
Inasimilarwayaswasdoneintheprevioussimulationexperimentwegeneratesequences{X ,...,X }
1 n
butin this case if {cp ,cp +l,...,cp ,cp +l} are the randomlyselected change-points,thenif
1 1 m/2 m/2
cp ≤i<cp +lforsomej ∈{1,...,m/2}thenX isdrawnfromaregionwithadeletionorasimple-
j j i
gaininoneallele(thisisa exclusive”or”),otherwiseisdrawnfromacopy-numberneutralsample. The
datasetchosenforoursimulationsislabeled”GSE29172”inacnr.
Tables4,5,6,7,8,9. Someofthegeneralobservationsmadeinthelastsubsectionremainstruehere.
InallstancesSaMe is thefasts methodbya largemarginandin termsof false-discoveryrateSaMe and
(cid:13)c 2010WILEY-VCHVerlagGmbH&Co.KGaA,Weinheim www.biometrical-journal.com
8 Pinheiroetal.anddd:Runningtitle
CBShaveasimilarperformance.PELThaveaverysmallnumberoffalse-discoverybutthisismainlydue
toitsinabilitytodetectsmallsegments.
Whenthealteredsegmentsweresamplefromthecopy-number1datasetallmethodsareabletodetect
morethan95%ofthepresentchange-points. Butwhenthealteredchange-pointsarefromadatasetwith
copy-number3 Table 7 shows that CBS has a slight better ability to detect change-pointsin relation to
SaMe,whilePELTdetectonlyasmallpercentageofthosechanges.InTables8and9weseethatwhenthe
alteredsegmentsareofsize50orgreaterallmethodsdetectmorethan95%ofthepresentchange-points.
Table4 ComparisonofSaMe,CBSandPELTwhenalteredsegmentsareofcopynumber1andoflength
25. SaMehasthefastestexecutiontimewhileCBShasaslightbetterchange-pointdetectionabilityand
PELTproducethesmallestnumberoffalse-discovery.
SaMe CBS PELT
t 0.0867 0.5887 0.4353
p10 0.9500 0.9900 0.9500
single
p5 0.9300 0.9750 0.9450
FP 0.2200 0.2300 0.0000
t 0.0769 0.6212 0.2840
p10 0.9875 1.0000 0.9275
double
p5 0.9775 0.9950 0.9200
FP 0.1300 0.1600 0.0000
t 0.0709 0.6473 0.2243
p10 0.9700 1.0000 0.9267
trip
p5 0.9517 0.9933 0.9200
FP 0.1800 0.1700 0.0000
Table5 ComparisonofSaMe,CBSandPELTwhenalteredsegmentsareofcopynumber1andoflength
50. SaMehasthefastestexecution,PELTproducethesmallestnumberoffalse-discoveryandallmethods
detectallpresentchange-points.
SaMe CBS PELT
t 0.0893 0.6596 0.4597
p10 1.0000 1.0000 1.0000
single
p5 0.9950 0.9950 0.9950
FP 0.2500 0.1600 0.0000
t 0.0750 0.6012 0.2763
p10 1.0000 1.0000 1.0000
double
p5 0.9925 0.9925 0.9950
FP 0.2200 0.1900 0.0000
t 0.0733 0.6312 0.2201
p10 1.0000 1.0000 1.0000
trip
p5 0.9917 0.9950 0.9933
FP 0.2200 0.2700 0.0000
(cid:13)c 2010WILEY-VCHVerlagGmbH&Co.KGaA,Weinheim www.biometrical-journal.com
BiometricalJournal52(2010)61 9
Table6 ComparisonofSaMe,CBSandPELTwhenalteredsegmentsareofcopynumber1andoflength
100.SimilarresultstothepreviousTable5.
SaMe CBS PELT
t 0.0870 0.5656 0.4160
p10 1.0000 1.0000 1.0000
single
p5 1.0000 0.9900 0.9850
FP 0.4000 0.1000 0.0000
t 0.0742 0.5435 0.2712
p10 1.0000 1.0000 1.0000
double
p5 0.9900 0.9875 0.9925
FP 0.2400 0.1800 0.0000
t 0.0752 0.6267 0.2169
p10 0.9967 1.0000 1.0000
trip
p5 0.9933 0.9917 0.9883
FP 0.2400 0.1700 0.0000
Table7 ComparisonofSaMe,CBSandPELTwhenalteredsegmentsareofcopynumber3andoflength
25.SaMehasthefastestexecutiontimewhileCBShasabetterchange-pointdetectionabilityspatiallyfor
thecasewhereonlyonealteredsegmentispresent. PELTproducethesmallestnumberoffalse-discovery
andwasparticularlyunabletodetectchange-points.
SaMe CBS PELT
t 0.0854 0.7368 0.4577
p10 0.7700 0.8450 0.3050
single
p5 0.7400 0.7700 0.3000
FP 0.3100 0.2100 0.0100
t 0.0852 1.0256 0.3019
p10 0.8800 0.8725 0.3450
double
p5 0.8550 0.8200 0.3325
FP 0.1400 0.1900 0.0400
t 0.0802 1.0088 0.2400
p10 0.8350 0.8850 0.3617
trip
p5 0.7867 0.8183 0.3517
FP 0.3300 0.3400 0.0300
5 Real Data Examples
In this section we show how the SaMe and other methods perform on real datasets. The first dataset
weconsider,availableathttp://penncnv.openbioinformatics.org,iscomprisedofafather-
mother-offspring trio produced with Illumina 550K platform. Following the analysis in NiuandZhang
(2012) we focusonthe offspringlog-R-ratiosequenceforchromosomes3, 11 and20. To each of these
sequencesweapplySaMe,CBSandPELT.TheparametersforSaMearek =25,k′ =20andpValuesof
(cid:13)c 2010WILEY-VCHVerlagGmbH&Co.KGaA,Weinheim www.biometrical-journal.com
10 Pinheiroetal.anddd:Runningtitle
Table8 ComparisonofSaMe,CBSandPELTwhenalteredsegmentsareofcopynumber3andoflength
50.SaMehasthefastestexecutiontimeandCBShasaslightlybetterchange-pointdetectionability.PELT
producethesmallestnumberoffalse-discovery.
SaMe CBS PELT
t 0.0862 0.6367 0.4112
p10 0.9550 0.9800 0.9550
single
p5 0.8800 0.8950 0.8550
FP 0.4700 0.2200 0.0300
t 0.0837 0.5872 0.2651
p10 0.9725 0.9950 0.9475
double
p5 0.9175 0.9200 0.9050
FP 0.3400 0.2300 0.1300
t 0.0731 0.6653 0.2263
p10 0.9733 0.9783 0.9717
trip
p5 0.9117 0.9200 0.9183
FP 0.3300 0.4000 0.0700
Table 9 Comparison of SaMe, CBS and PELT when altered segments are of copy number 3 and of
length100. SaMehasthefastestexecutiontime. Withrespectofothercriteriathemethodshaveasimilar
behavior.
SaMe CBS PELT
t 0.0870 0.6285 0.3817
p10 0.9800 0.9900 0.9750
single
p5 0.9050 0.9300 0.9000
FP 0.2300 0.1300 0.0500
t 0.0797 0.6002 0.2822
p10 0.9775 0.9850 0.9875
double
p5 0.9275 0.9275 0.9225
FP 0.2600 0.2100 0.0500
t 0.0728 0.6383 0.2265
p10 0.9817 0.9833 0.9883
trip
p5 0.9217 0.9217 0.9350
FP 0.3100 0.3300 0.0700
0.01foranytestperformed. We chooseCBSparametersinthesamewaywedidintheprevioussection
with the exception that here we adopt a change-point pruning procedure that removes a change-point.
PELTparametersarethesameonesasintheprevioussectionandweconsiderseparatelytheBICandAIC
criterions.
Table 10 shows the number of change-points detected by each method. Our results are very similar
to those found in NiuandZhang (2012) where the execution of PennCNV Wangetal. (2007) found 2,
(cid:13)c 2010WILEY-VCHVerlagGmbH&Co.KGaA,Weinheim www.biometrical-journal.com