R Under development (unstable) (2024-02-16 r85931 ucrt) -- "Unsuffered Consequences" Copyright (C) 2024 The R Foundation for Statistical Computing Platform: x86_64-w64-mingw32/x64 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > library("PSCBS") PSCBS v0.67.0 successfully loaded. See ?PSCBS for help. > > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > # Load SNP microarray data > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > data <- PSCBS::exampleData("paired.chr01") > str(data) 'data.frame': 73346 obs. of 6 variables: $ chromosome: int 1 1 1 1 1 1 1 1 1 1 ... $ x : int 1145994 2224111 2319424 2543484 2926730 2941694 3084986 3155127 3292731 3695086 ... $ CT : num 1.625 1.071 1.406 1.18 0.856 ... $ betaT : num 0.757 0.771 0.834 0.778 0.229 ... $ CN : num 2.36 2.13 2.59 1.93 1.71 ... $ betaN : num 0.827 0.875 0.887 0.884 0.103 ... > > > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > # Paired PSCBS segmentation > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > # Drop single-locus outliers > dataS <- dropSegmentationOutliers(data) > > # Find centromere > gaps <- findLargeGaps(dataS, minLength=2e6) > knownSegments <- gapsToSegments(gaps) > > > # Run light-weight tests by default > if (Sys.getenv("_R_CHECK_FULL_") == "") { + # Use only every 5th data point + dataS <- dataS[seq(from=1, to=nrow(data), by=5),] + # Number of segments (for assertion) + nSegs <- 4L + # Number of bootstrap samples (see below) + B <- 100L + } else { + # Full tests + nSegs <- 11L + B <- 1000L + } > > str(dataS) 'data.frame': 14670 obs. of 6 variables: $ chromosome: int 1 1 1 1 1 1 1 1 1 1 ... $ x : int 1145994 2941694 3710825 4240737 4276892 4464544 4714611 5095111 5034491 5158416 ... $ CT : num 1.63 1.4 1.41 1.17 1.16 ... $ betaT : num 0.7574 0.169 0.2357 0.2604 0.0576 ... $ CN : num 2.36 2.13 2.26 2.01 2.32 ... $ betaN : num 0.8274 0.5072 0.1671 0.1609 0.0421 ... > > # Paired PSCBS segmentation > fit <- segmentByPairedPSCBS(dataS, knownSegments=knownSegments, + seed=0xBEEF, verbose=-10) Segmenting paired tumor-normal signals using Paired PSCBS... Calling genotypes from normal allele B fractions... num [1:14670] 0.8274 0.5072 0.1671 0.1609 0.0421 ... Called genotypes: num [1:14670] 1 0.5 0 0 0 0 1 0 1 0.5 ... - attr(*, "modelFit")=List of 1 ..$ :List of 7 .. ..$ flavor : chr "density" .. ..$ cn : int 2 .. ..$ nbrOfGenotypeGroups: int 3 .. ..$ tau : num [1:2] 0.315 0.677 .. ..$ n : int 14640 .. ..$ fit :Classes 'PeaksAndValleys' and 'data.frame': 5 obs. of 3 variables: .. .. ..$ type : chr [1:5] "peak" "valley" "peak" "valley" ... .. .. ..$ x : num [1:5] 0.104 0.315 0.499 0.677 0.885 .. .. ..$ density: num [1:5] 1.479 0.522 1.056 0.551 1.453 .. ..$ fitValleys :Classes 'PeaksAndValleys' and 'data.frame': 2 obs. of 3 variables: .. .. ..$ type : chr [1:2] "valley" "valley" .. .. ..$ x : num [1:2] 0.315 0.677 .. .. ..$ density: num [1:2] 0.522 0.551 ..- attr(*, "class")= chr [1:2] "NaiveGenotypeModelFit" "list" muN 0 0.5 1 5221 4198 5251 Calling genotypes from normal allele B fractions...done Normalizing betaT using betaN (TumorBoost)... Normalized BAFs: num [1:14670] 0.9301 0.1667 0.0685 0.0995 0.0155 ... - attr(*, "modelFit")=List of 5 ..$ method : chr "normalizeTumorBoost" ..$ flavor : chr "v4" ..$ delta : num [1:14670] -0.17264 0.00239 0.1671 0.16085 0.04213 ... .. ..- attr(*, "modelFit")=List of 1 .. .. ..$ :List of 7 .. .. .. ..$ flavor : chr "density" .. .. .. ..$ cn : int 2 .. .. .. ..$ nbrOfGenotypeGroups: int 3 .. .. .. ..$ tau : num [1:2] 0.315 0.677 .. .. .. ..$ n : int 14640 .. .. .. ..$ fit :Classes 'PeaksAndValleys' and 'data.frame': 5 obs. of 3 variables: .. .. .. .. ..$ type : chr [1:5] "peak" "valley" "peak" "valley" ... .. .. .. .. ..$ x : num [1:5] 0.104 0.315 0.499 0.677 0.885 .. .. .. .. ..$ density: num [1:5] 1.479 0.522 1.056 0.551 1.453 .. .. .. ..$ fitValleys :Classes 'PeaksAndValleys' and 'data.frame': 2 obs. of 3 variables: .. .. .. .. ..$ type : chr [1:2] "valley" "valley" .. .. .. .. ..$ x : num [1:2] 0.315 0.677 .. .. .. .. ..$ density: num [1:2] 0.522 0.551 .. .. ..- attr(*, "class")= chr [1:2] "NaiveGenotypeModelFit" "list" ..$ preserveScale: logi FALSE ..$ scaleFactor : num NA Normalizing betaT using betaN (TumorBoost)...done Setup up data... 'data.frame': 14670 obs. of 7 variables: $ chromosome: int 1 1 1 1 1 1 1 1 1 1 ... $ x : num 1145994 2941694 3710825 4240737 4276892 ... $ CT : num 1.63 1.4 1.41 1.17 1.16 ... $ betaT : num 0.7574 0.169 0.2357 0.2604 0.0576 ... $ betaTN : num 0.9301 0.1667 0.0685 0.0995 0.0155 ... ..- attr(*, "modelFit")=List of 5 .. ..$ method : chr "normalizeTumorBoost" .. ..$ flavor : chr "v4" .. ..$ delta : num [1:14670] -0.17264 0.00239 0.1671 0.16085 0.04213 ... .. .. ..- attr(*, "modelFit")=List of 1 .. .. .. ..$ :List of 7 .. .. .. .. ..$ flavor : chr "density" .. .. .. .. ..$ cn : int 2 .. .. .. .. ..$ nbrOfGenotypeGroups: int 3 .. .. .. .. ..$ tau : num [1:2] 0.315 0.677 .. .. .. .. ..$ n : int 14640 .. .. .. .. ..$ fit :Classes 'PeaksAndValleys' and 'data.frame': 5 obs. of 3 variables: .. .. .. .. .. ..$ type : chr [1:5] "peak" "valley" "peak" "valley" ... .. .. .. .. .. ..$ x : num [1:5] 0.104 0.315 0.499 0.677 0.885 .. .. .. .. .. ..$ density: num [1:5] 1.479 0.522 1.056 0.551 1.453 .. .. .. .. ..$ fitValleys :Classes 'PeaksAndValleys' and 'data.frame': 2 obs. of 3 variables: .. .. .. .. .. ..$ type : chr [1:2] "valley" "valley" .. .. .. .. .. ..$ x : num [1:2] 0.315 0.677 .. .. .. .. .. ..$ density: num [1:2] 0.522 0.551 .. .. .. ..- attr(*, "class")= chr [1:2] "NaiveGenotypeModelFit" "list" .. ..$ preserveScale: logi FALSE .. ..$ scaleFactor : num NA $ betaN : num 0.8274 0.5072 0.1671 0.1609 0.0421 ... $ muN : num 1 0.5 0 0 0 0 1 0 1 0.5 ... ..- attr(*, "modelFit")=List of 1 .. ..$ :List of 7 .. .. ..$ flavor : chr "density" .. .. ..$ cn : int 2 .. .. ..$ nbrOfGenotypeGroups: int 3 .. .. ..$ tau : num [1:2] 0.315 0.677 .. .. ..$ n : int 14640 .. .. ..$ fit :Classes 'PeaksAndValleys' and 'data.frame': 5 obs. of 3 variables: .. .. .. ..$ type : chr [1:5] "peak" "valley" "peak" "valley" ... .. .. .. ..$ x : num [1:5] 0.104 0.315 0.499 0.677 0.885 .. .. .. ..$ density: num [1:5] 1.479 0.522 1.056 0.551 1.453 .. .. ..$ fitValleys :Classes 'PeaksAndValleys' and 'data.frame': 2 obs. of 3 variables: .. .. .. ..$ type : chr [1:2] "valley" "valley" .. .. .. ..$ x : num [1:2] 0.315 0.677 .. .. .. ..$ density: num [1:2] 0.522 0.551 .. ..- attr(*, "class")= chr [1:2] "NaiveGenotypeModelFit" "list" Setup up data...done Dropping loci for which TCNs are missing... Number of loci dropped: 12 Dropping loci for which TCNs are missing...done Ordering data along genome... 'data.frame': 14658 obs. of 7 variables: $ chromosome: int 1 1 1 1 1 1 1 1 1 1 ... $ x : num 554484 730720 782343 878522 916294 ... $ CT : num 1.88 1.8 1.59 1.64 1.53 ... $ betaT : num 0.0646 0.1672 0.9284 0.113 0.7209 ... $ betaTN : num -0.0515 -0.1172 1.0194 0.031 0.8604 ... $ betaN : num 0.116 0.284 0.909 0.082 0.86 ... $ muN : num 0 0 1 0 1 1 1 0 1 0.5 ... Ordering data along genome...done Keeping only current chromosome for 'knownSegments'... Chromosome: 1 Known segments for this chromosome: chromosome start end length 1 1 -Inf 120992603 Inf 2 1 120992604 141510002 20517398 3 1 141510003 Inf Inf Keeping only current chromosome for 'knownSegments'...done alphaTCN: 0.009 alphaDH: 0.001 Number of loci: 14658 Calculating DHs... Number of SNPs: 14658 Number of heterozygous SNPs: 4196 (28.63%) Normalized DHs: num [1:14658] NA NA NA NA NA ... Calculating DHs...done Random seed temporarily set (seed=c(48879), kind="L'Ecuyer-CMRG") Produced 2 seeds from this stream for future usage Identification of change points by total copy numbers... Segmenting by CBS... Chromosome: 1 Segmenting multiple segments on current chromosome... Number of segments: 3 Random seed temporarily set (seed=c(10407, 1066287653, -51199871, 161854402, -1995183193, 1503453565, -747102133), kind="L'Ecuyer-CMRG") Produced 3 seeds from this stream for future usage Segmenting by CBS... Chromosome: 1 Random seed temporarily set (seed=c(10407, 1797822437, 388243314, 91406689, -519578635, -1381924756, 1089253019), kind="L'Ecuyer-CMRG") Segmenting by CBS...done Segmenting by CBS... Chromosome: 1 Random seed temporarily set (seed=c(10407, -1924040949, -1632234809, -437763632, -1464377300, 676654412, 2080370711), kind="L'Ecuyer-CMRG") Segmenting by CBS...done Segmenting multiple segments on current chromosome...done Segmenting by CBS...done List of 4 $ data :'data.frame': 14658 obs. of 4 variables: ..$ chromosome: int [1:14658] 1 1 1 1 1 1 1 1 1 1 ... ..$ x : num [1:14658] 554484 730720 782343 878522 916294 ... ..$ y : num [1:14658] 1.88 1.8 1.59 1.64 1.53 ... ..$ index : int [1:14658] 1 2 3 4 5 6 7 8 9 10 ... $ output :'data.frame': 4 obs. of 6 variables: ..$ sampleName: chr [1:4] NA NA NA NA ..$ chromosome: int [1:4] 1 1 1 1 ..$ start : num [1:4] 5.54e+05 1.21e+08 1.42e+08 1.85e+08 ..$ end : num [1:4] 1.21e+08 1.42e+08 1.85e+08 2.47e+08 ..$ nbrOfLoci : int [1:4] 7586 0 2681 4391 ..$ mean : num [1:4] 1.39 NA 2.07 2.63 $ segRows:'data.frame': 4 obs. of 2 variables: ..$ startRow: int [1:4] 1 NA 7587 10268 ..$ endRow : int [1:4] 7586 NA 10267 14658 $ params :List of 5 ..$ alpha : num 0.009 ..$ undo : num 0 ..$ joinSegments : logi TRUE ..$ knownSegments:'data.frame': 4 obs. of 3 variables: .. ..$ chromosome: int [1:4] 1 1 2 1 .. ..$ start : num [1:4] -Inf -Inf -Inf 1.42e+08 .. ..$ end : num [1:4] 1.21e+08 Inf Inf Inf ..$ seed : int [1:7] 10407 1066287653 -51199871 161854402 -1995183193 1503453565 -747102133 - attr(*, "class")= chr [1:2] "CBS" "AbstractCBS" - attr(*, "processingTime")= 'proc_time' Named num [1:5] 0.23 0 0.23 NA NA ..- attr(*, "names")= chr [1:5] "user.self" "sys.self" "elapsed" "user.child" ... - attr(*, "pkgDetails")= chr "DNAcopy v1.76.0" - attr(*, "randomSeed")= int [1:7] 10407 1797822437 388243314 91406689 -519578635 -1381924756 1089253019 Identification of change points by total copy numbers...done Restructure TCN segmentation results... chromosome tcnStart tcnEnd tcnNbrOfLoci tcnMean 1 1 554484 120992603 7586 1.3853 2 1 120992604 141510002 0 NA 3 1 141510003 185449813 2681 2.0689 4 1 185449813 247137334 4391 2.6341 Number of TCN segments: 4 Restructure TCN segmentation results...done Total CN segment #1 ([ 554484,1.20993e+08]) of 4... Number of TCN loci in segment: 7586 Locus data for TCN segment: 'data.frame': 7586 obs. of 9 variables: $ chromosome: int 1 1 1 1 1 1 1 1 1 1 ... $ x : num 554484 730720 782343 878522 916294 ... $ CT : num 1.88 1.8 1.59 1.64 1.53 ... $ betaT : num 0.0646 0.1672 0.9284 0.113 0.7209 ... $ betaTN : num -0.0515 -0.1172 1.0194 0.031 0.8604 ... $ betaN : num 0.116 0.284 0.909 0.082 0.86 ... $ muN : num 0 0 1 0 1 1 1 0 1 0.5 ... $ index : int 1 2 3 4 5 6 7 8 9 10 ... $ rho : num NA NA NA NA NA ... Number of loci: 7586 Number of SNPs: 2108 (27.79%) Number of heterozygous SNPs: 2108 (100.00%) Chromosome: 1 Segmenting DH signals... Segmenting by CBS... Chromosome: 1 Random seed temporarily set (seed=c(10407, 1797822437, 388243314, 91406689, -519578635, -1381924756, 1089253019), kind="L'Ecuyer-CMRG") Segmenting by CBS...done List of 4 $ data :'data.frame': 7586 obs. of 4 variables: ..$ chromosome: int [1:7586] 1 1 1 1 1 1 1 1 1 1 ... ..$ x : num [1:7586] 554484 730720 782343 878522 916294 ... ..$ y : num [1:7586] NA NA NA NA NA ... ..$ index : int [1:7586] 1 2 3 4 5 6 7 8 9 10 ... $ output :'data.frame': 1 obs. of 6 variables: ..$ sampleName: chr NA ..$ chromosome: int 1 ..$ start : num 554484 ..$ end : num 1.21e+08 ..$ nbrOfLoci : int 2108 ..$ mean : num 0.512 $ segRows:'data.frame': 1 obs. of 2 variables: ..$ startRow: int 10 ..$ endRow : int 7574 $ params :List of 5 ..$ alpha : num 0.001 ..$ undo : num 0 ..$ joinSegments : logi TRUE ..$ knownSegments:'data.frame': 1 obs. of 3 variables: .. ..$ chromosome: int 1 .. ..$ start : num 554484 .. ..$ end : num 1.21e+08 ..$ seed : int [1:7] 10407 1797822437 388243314 91406689 -519578635 -1381924756 1089253019 - attr(*, "class")= chr [1:2] "CBS" "AbstractCBS" - attr(*, "processingTime")= 'proc_time' Named num [1:5] 0.1 0 0.1 NA NA ..- attr(*, "names")= chr [1:5] "user.self" "sys.self" "elapsed" "user.child" ... - attr(*, "pkgDetails")= chr "DNAcopy v1.76.0" - attr(*, "randomSeed")= int [1:7] 10407 1797822437 388243314 91406689 -519578635 -1381924756 1089253019 DH segmentation (locally-indexed) rows: startRow endRow 1 10 7574 int [1:7586] 1 2 3 4 5 6 7 8 9 10 ... DH segmentation rows: startRow endRow 1 10 7574 Segmenting DH signals...done DH segmentation table: dhStart dhEnd dhNbrOfLoci dhMean 1 554484 120992603 2108 0.5116 startRow endRow 1 10 7574 Rows: [1] 1 TCN segmentation rows: startRow endRow 1 1 7586 TCN and DH segmentation rows: startRow endRow 1 1 7586 startRow endRow 1 10 7574 NULL TCN segmentation (expanded) rows: startRow endRow 1 1 7586 TCN and DH segmentation rows: startRow endRow 1 1 7586 2 NA NA 3 7587 10267 4 10268 14658 startRow endRow 1 10 7574 startRow endRow 1 1 7586 Total CN segmentation table (expanded): chromosome tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs tcnNbrOfHets 1 1 554484 120992603 7586 1.3853 2108 2108 (TCN,DH) segmentation for one total CN segment: tcnId dhId chromosome tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs 1 1 1 1 554484 120992603 7586 1.3853 2108 tcnNbrOfHets dhStart dhEnd dhNbrOfLoci dhMean 1 2108 554484 120992603 2108 0.5116 Total CN segment #1 ([ 554484,1.20993e+08]) of 4...done Total CN segment #2 ([1.20993e+08,1.4151e+08]) of 4... Number of TCN loci in segment: 0 Locus data for TCN segment: 'data.frame': 0 obs. of 9 variables: $ chromosome: int $ x : num $ CT : num $ betaT : num $ betaTN : num $ betaN : num $ muN : num $ index : int $ rho : num Number of loci: 0 Number of SNPs: 0 (NaN%) Number of heterozygous SNPs: 0 (NaN%) Chromosome: 1 Segmenting DH signals... Segmenting by CBS... Chromosome: NA Random seed temporarily set (seed=c(10407, 1797822437, 388243314, 91406689, -519578635, -1381924756, 1089253019), kind="L'Ecuyer-CMRG") Segmenting by CBS...done List of 4 $ data :'data.frame': 0 obs. of 4 variables: ..$ chromosome: int(0) ..$ x : num(0) ..$ y : num(0) ..$ index : int(0) $ output :'data.frame': 0 obs. of 6 variables: ..$ sampleName: chr(0) ..$ chromosome: num(0) ..$ start : num(0) ..$ end : num(0) ..$ nbrOfLoci : int(0) ..$ mean : num(0) $ segRows:'data.frame': 0 obs. of 2 variables: ..$ startRow: int(0) ..$ endRow : int(0) $ params :List of 5 ..$ alpha : num 0.001 ..$ undo : num 0 ..$ joinSegments : logi TRUE ..$ knownSegments:'data.frame': 0 obs. of 3 variables: .. ..$ chromosome: int(0) .. ..$ start : num(0) .. ..$ end : num(0) ..$ seed : int [1:7] 10407 1797822437 388243314 91406689 -519578635 -1381924756 1089253019 - attr(*, "class")= chr [1:2] "CBS" "AbstractCBS" - attr(*, "processingTime")= 'proc_time' Named num [1:5] 0 0 0 NA NA ..- attr(*, "names")= chr [1:5] "user.self" "sys.self" "elapsed" "user.child" ... - attr(*, "pkgDetails")= chr "DNAcopy v1.76.0" - attr(*, "randomSeed")= int [1:7] 10407 1797822437 388243314 91406689 -519578635 -1381924756 1089253019 DH segmentation (locally-indexed) rows: [1] startRow endRow <0 rows> (or 0-length row.names) int(0) DH segmentation rows: [1] startRow endRow <0 rows> (or 0-length row.names) Segmenting DH signals...done DH segmentation table: dhStart dhEnd dhNbrOfLoci dhMean NA NA NA NA NA startRow endRow NA NA NA Rows: [1] 2 TCN segmentation rows: startRow endRow 2 NA NA TCN and DH segmentation rows: startRow endRow 2 NA NA startRow endRow NA NA NA startRow endRow 1 1 7586 TCN segmentation (expanded) rows: startRow endRow 1 1 7586 2 NA NA TCN and DH segmentation rows: startRow endRow 1 1 7586 2 NA NA 3 7587 10267 4 10268 14658 startRow endRow 1 10 7574 2 NA NA startRow endRow 1 1 7586 2 NA NA Total CN segmentation table (expanded): chromosome tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs tcnNbrOfHets 2 1 120992604 141510002 0 NA 0 0 (TCN,DH) segmentation for one total CN segment: tcnId dhId chromosome tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs 2 2 1 1 120992604 141510002 0 NA 0 tcnNbrOfHets dhStart dhEnd dhNbrOfLoci dhMean 2 0 NA NA NA NA Total CN segment #2 ([1.20993e+08,1.4151e+08]) of 4...done Total CN segment #3 ([1.4151e+08,1.8545e+08]) of 4... Number of TCN loci in segment: 2681 Locus data for TCN segment: 'data.frame': 2681 obs. of 9 variables: $ chromosome: int 1 1 1 1 1 1 1 1 1 1 ... $ x : num 1.43e+08 1.43e+08 1.43e+08 1.43e+08 1.44e+08 ... $ CT : num 2.27 1.55 1.47 1.5 1.81 ... $ betaT : num 0.34 0.55 0.048 0.813 0.575 ... $ betaTN : num 0.441 0.629 -0.05 0.958 0.872 ... $ betaN : num 0.3851 0.3933 0.0981 0.8552 0.7028 ... $ muN : num 0.5 0.5 0 1 1 1 1 0.5 1 1 ... $ index : int 7587 7588 7589 7590 7591 7592 7593 7594 7595 7596 ... $ rho : num 0.117 0.258 NA NA NA ... Number of loci: 2681 Number of SNPs: 777 (28.98%) Number of heterozygous SNPs: 777 (100.00%) Chromosome: 1 Segmenting DH signals... Segmenting by CBS... Chromosome: 1 Random seed temporarily set (seed=c(10407, 1797822437, 388243314, 91406689, -519578635, -1381924756, 1089253019), kind="L'Ecuyer-CMRG") Segmenting by CBS...done List of 4 $ data :'data.frame': 2681 obs. of 4 variables: ..$ chromosome: int [1:2681] 1 1 1 1 1 1 1 1 1 1 ... ..$ x : num [1:2681] 1.43e+08 1.43e+08 1.43e+08 1.43e+08 1.44e+08 ... ..$ y : num [1:2681] 0.117 0.258 NA NA NA ... ..$ index : int [1:2681] 1 2 3 4 5 6 7 8 9 10 ... $ output :'data.frame': 1 obs. of 6 variables: ..$ sampleName: chr NA ..$ chromosome: int 1 ..$ start : num 1.42e+08 ..$ end : num 1.85e+08 ..$ nbrOfLoci : int 777 ..$ mean : num 0.0973 $ segRows:'data.frame': 1 obs. of 2 variables: ..$ startRow: int 1 ..$ endRow : int 2677 $ params :List of 5 ..$ alpha : num 0.001 ..$ undo : num 0 ..$ joinSegments : logi TRUE ..$ knownSegments:'data.frame': 1 obs. of 3 variables: .. ..$ chromosome: int 1 .. ..$ start : num 1.42e+08 .. ..$ end : num 1.85e+08 ..$ seed : int [1:7] 10407 1797822437 388243314 91406689 -519578635 -1381924756 1089253019 - attr(*, "class")= chr [1:2] "CBS" "AbstractCBS" - attr(*, "processingTime")= 'proc_time' Named num [1:5] 0.02 0 0.02 NA NA ..- attr(*, "names")= chr [1:5] "user.self" "sys.self" "elapsed" "user.child" ... - attr(*, "pkgDetails")= chr "DNAcopy v1.76.0" - attr(*, "randomSeed")= int [1:7] 10407 1797822437 388243314 91406689 -519578635 -1381924756 1089253019 DH segmentation (locally-indexed) rows: startRow endRow 1 1 2677 int [1:2681] 7587 7588 7589 7590 7591 7592 7593 7594 7595 7596 ... DH segmentation rows: startRow endRow 1 7587 10263 Segmenting DH signals...done DH segmentation table: dhStart dhEnd dhNbrOfLoci dhMean 1 141510003 185449813 777 0.0973 startRow endRow 1 7587 10263 Rows: [1] 3 TCN segmentation rows: startRow endRow 3 7587 10267 TCN and DH segmentation rows: startRow endRow 3 7587 10267 startRow endRow 1 7587 10263 startRow endRow 1 1 7586 2 NA NA TCN segmentation (expanded) rows: startRow endRow 1 1 7586 2 NA NA 3 7587 10267 TCN and DH segmentation rows: startRow endRow 1 1 7586 2 NA NA 3 7587 10267 4 10268 14658 startRow endRow 1 10 7574 2 NA NA 3 7587 10263 startRow endRow 1 1 7586 2 NA NA 3 7587 10267 Total CN segmentation table (expanded): chromosome tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs tcnNbrOfHets 3 1 141510003 185449813 2681 2.0689 777 777 (TCN,DH) segmentation for one total CN segment: tcnId dhId chromosome tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs 3 3 1 1 141510003 185449813 2681 2.0689 777 tcnNbrOfHets dhStart dhEnd dhNbrOfLoci dhMean 3 777 141510003 185449813 777 0.0973 Total CN segment #3 ([1.4151e+08,1.8545e+08]) of 4...done Total CN segment #4 ([1.8545e+08,2.47137e+08]) of 4... Number of TCN loci in segment: 4391 Locus data for TCN segment: 'data.frame': 4391 obs. of 9 variables: $ chromosome: int 1 1 1 1 1 1 1 1 1 1 ... $ x : num 1.85e+08 1.85e+08 1.85e+08 1.86e+08 1.86e+08 ... $ CT : num 2.93 2.15 2.82 2.93 2.46 ... $ betaT : num 0.0811 0.5154 0.9473 0.3734 0.7506 ... $ betaTN : num -0.169 0.609 1.028 0.525 0.968 ... $ betaN : num 0.25 0.38 0.919 0.34 0.783 ... $ muN : num 0 0.5 1 0.5 1 1 0 1 0 1 ... $ index : int 10268 10269 10270 10271 10272 10273 10274 10275 10276 10277 ... $ rho : num NA 0.2186 NA 0.0503 NA ... Number of loci: 4391 Number of SNPs: 1311 (29.86%) Number of heterozygous SNPs: 1311 (100.00%) Chromosome: 1 Segmenting DH signals... Segmenting by CBS... Chromosome: 1 Random seed temporarily set (seed=c(10407, 1797822437, 388243314, 91406689, -519578635, -1381924756, 1089253019), kind="L'Ecuyer-CMRG") Segmenting by CBS...done List of 4 $ data :'data.frame': 4391 obs. of 4 variables: ..$ chromosome: int [1:4391] 1 1 1 1 1 1 1 1 1 1 ... ..$ x : num [1:4391] 1.85e+08 1.85e+08 1.85e+08 1.86e+08 1.86e+08 ... ..$ y : num [1:4391] NA 0.2186 NA 0.0503 NA ... ..$ index : int [1:4391] 1 2 3 4 5 6 7 8 9 10 ... $ output :'data.frame': 1 obs. of 6 variables: ..$ sampleName: chr NA ..$ chromosome: int 1 ..$ start : num 1.85e+08 ..$ end : num 2.47e+08 ..$ nbrOfLoci : int 1311 ..$ mean : num 0.23 $ segRows:'data.frame': 1 obs. of 2 variables: ..$ startRow: int 2 ..$ endRow : int 4388 $ params :List of 5 ..$ alpha : num 0.001 ..$ undo : num 0 ..$ joinSegments : logi TRUE ..$ knownSegments:'data.frame': 1 obs. of 3 variables: .. ..$ chromosome: int 1 .. ..$ start : num 1.85e+08 .. ..$ end : num 2.47e+08 ..$ seed : int [1:7] 10407 1797822437 388243314 91406689 -519578635 -1381924756 1089253019 - attr(*, "class")= chr [1:2] "CBS" "AbstractCBS" - attr(*, "processingTime")= 'proc_time' Named num [1:5] 0.03 0 0.03 NA NA ..- attr(*, "names")= chr [1:5] "user.self" "sys.self" "elapsed" "user.child" ... - attr(*, "pkgDetails")= chr "DNAcopy v1.76.0" - attr(*, "randomSeed")= int [1:7] 10407 1797822437 388243314 91406689 -519578635 -1381924756 1089253019 DH segmentation (locally-indexed) rows: startRow endRow 1 2 4388 int [1:4391] 10268 10269 10270 10271 10272 10273 10274 10275 10276 10277 ... DH segmentation rows: startRow endRow 1 10269 14655 Segmenting DH signals...done DH segmentation table: dhStart dhEnd dhNbrOfLoci dhMean 1 185449813 247137334 1311 0.2295 startRow endRow 1 10269 14655 Rows: [1] 4 TCN segmentation rows: startRow endRow 4 10268 14658 TCN and DH segmentation rows: startRow endRow 4 10268 14658 startRow endRow 1 10269 14655 startRow endRow 1 1 7586 2 NA NA 3 7587 10267 TCN segmentation (expanded) rows: startRow endRow 1 1 7586 2 NA NA 3 7587 10267 4 10268 14658 TCN and DH segmentation rows: startRow endRow 1 1 7586 2 NA NA 3 7587 10267 4 10268 14658 startRow endRow 1 10 7574 2 NA NA 3 7587 10263 4 10269 14655 startRow endRow 1 1 7586 2 NA NA 3 7587 10267 4 10268 14658 Total CN segmentation table (expanded): chromosome tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs tcnNbrOfHets 4 1 185449813 247137334 4391 2.6341 1311 1311 (TCN,DH) segmentation for one total CN segment: tcnId dhId chromosome tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs 4 4 1 1 185449813 247137334 4391 2.6341 1311 tcnNbrOfHets dhStart dhEnd dhNbrOfLoci dhMean 4 1311 185449813 247137334 1311 0.2295 Total CN segment #4 ([1.8545e+08,2.47137e+08]) of 4...done chromosome tcnId dhId tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs 1 1 1 1 554484 120992603 7586 1.3853 2108 2 1 2 1 120992604 141510002 0 NA 0 3 1 3 1 141510003 185449813 2681 2.0689 777 4 1 4 1 185449813 247137334 4391 2.6341 1311 tcnNbrOfHets dhStart dhEnd dhNbrOfLoci dhMean 1 2108 554484 120992603 2108 0.5116 2 0 NA NA NA NA 3 777 141510003 185449813 777 0.0973 4 1311 185449813 247137334 1311 0.2295 Calculating (C1,C2) per segment... Calculating (C1,C2) per segment...done Number of segments: 4 Segmenting paired tumor-normal signals using Paired PSCBS...done Post-segmenting TCNs... Number of segments: 4 Number of chromosomes: 1 [1] 1 Chromosome 1 ('chr01') of 1... Rows: [1] 1 2 3 4 Number of segments: 4 TCN segment #1 ('1') of 4... Nothing todo. Only one DH segmentation. Skipping. TCN segment #1 ('1') of 4...done TCN segment #2 ('2') of 4... Nothing todo. Only one DH segmentation. Skipping. TCN segment #2 ('2') of 4...done TCN segment #3 ('3') of 4... Nothing todo. Only one DH segmentation. Skipping. TCN segment #3 ('3') of 4...done TCN segment #4 ('4') of 4... Nothing todo. Only one DH segmentation. Skipping. TCN segment #4 ('4') of 4...done Chromosome 1 ('chr01') of 1...done Update (C1,C2) per segment... Update (C1,C2) per segment...done Post-segmenting TCNs...done chromosome tcnId dhId tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs 1 1 1 1 554484 120992603 7586 1.3853 2108 2 1 2 1 120992604 141510002 0 NA 0 3 1 3 1 141510003 185449813 2681 2.0689 777 4 1 4 1 185449813 247137334 4391 2.6341 1311 tcnNbrOfHets dhStart dhEnd dhNbrOfLoci dhMean c1Mean c2Mean 1 2108 554484 120992603 2108 0.5116 0.3382903 1.047010 2 0 NA NA NA NA NA NA 3 777 141510003 185449813 777 0.0973 0.9337980 1.135102 4 1311 185449813 247137334 1311 0.2295 1.0147870 1.619313 chromosome tcnId dhId tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs 1 1 1 1 554484 120992603 7586 1.3853 2108 2 1 2 1 120992604 141510002 0 NA 0 3 1 3 1 141510003 185449813 2681 2.0689 777 4 1 4 1 185449813 247137334 4391 2.6341 1311 tcnNbrOfHets dhStart dhEnd dhNbrOfLoci dhMean c1Mean c2Mean 1 2108 554484 120992603 2108 0.5116 0.3382903 1.047010 2 0 NA NA NA NA NA NA 3 777 141510003 185449813 777 0.0973 0.9337980 1.135102 4 1311 185449813 247137334 1311 0.2295 1.0147870 1.619313 > print(fit) chromosome tcnId dhId start end tcnNbrOfLoci tcnMean tcnNbrOfSNPs 1 1 1 1 554484 120992603 7586 1.3853 2108 2 1 2 1 120992604 141510002 0 NA 0 3 1 3 1 141510003 185449813 2681 2.0689 777 4 1 4 1 185449813 247137334 4391 2.6341 1311 tcnNbrOfHets dhNbrOfLoci dhMean c1Mean c2Mean 1 2108 2108 0.5116 0.3382903 1.047010 2 0 NA NA NA NA 3 777 777 0.0973 0.9337980 1.135102 4 1311 1311 0.2295 1.0147870 1.619313 > > # Plot results > plotTracks(fit) > > # Sanity check > stopifnot(nbrOfSegments(fit) == nSegs) > > > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > # Bootstrap segment level estimates > # (used by the AB caller, which, if skipped here, > # will do it automatically) > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > fit <- bootstrapTCNandDHByRegion(fit, B=B, verbose=-10) Resample (TCN,DH) signals and re-estimate summaries for segment & changepoint... Already done? tcn_2.5% tcn_5% tcn_95% tcn_97.5% dh_2.5% dh_5% dh_95% dh_97.5% FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE c1_2.5% c1_5% c1_95% c1_97.5% c2_2.5% c2_5% c2_95% c2_97.5% FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE Bootstrapping (TCN,DH,C1,C2) segment mean levels... Identifying heterozygous & homozygous SNPs and non-polymorphic loci... Number of loci: 14658 Number of SNPs: 4196 Number of non-SNPs: 10462 Identifying heterozygous & homozygous SNPs and non-polymorphic loci...done num [1:4, 1:100, 1:4] NA NA NA NA NA NA NA NA NA NA ... - attr(*, "dimnames")=List of 3 ..$ : NULL ..$ : NULL ..$ : chr [1:4] "tcn" "dh" "c1" "c2" Segment #1 (chr 1, tcnId=1, dhId=1) of 4... chromosome tcnId dhId tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs 1 1 1 1 554484 120992603 7586 1.3853 2108 tcnNbrOfHets dhStart dhEnd dhNbrOfLoci dhMean c1Mean c2Mean 1 2108 554484 120992603 2108 0.5116 0.3382903 1.04701 Number of TCNs: 7586 Number of DHs: 2108 int [1:7586] 1 2 3 4 5 6 7 8 9 10 ... int [1:14658] 1 2 3 4 5 6 7 8 9 10 ... int [1:7586] 1 2 3 4 5 6 7 8 9 10 ... Identify loci used to bootstrap DH means... Heterozygous SNPs to resample for DH: int [1:2108] 10 12 24 28 31 33 34 39 46 48 ... Identify loci used to bootstrap DH means...done Identify loci used to bootstrap TCN means... SNPs: int [1:2108] 10 12 24 28 31 33 34 39 46 48 ... Non-polymorphic loci: int [1:5478] 1 2 3 4 5 6 7 8 9 11 ... Heterozygous SNPs to resample for TCN: int [1:2108] 10 12 24 28 31 33 34 39 46 48 ... Homozygous SNPs to resample for TCN: int(0) Non-polymorphic loci to resample for TCN: int [1:5478] 1 2 3 4 5 6 7 8 9 11 ... Heterozygous SNPs with non-DH to resample for TCN: int(0) Loci to resample for TCN: int [1:7586] 1 2 3 4 5 6 7 8 9 10 ... Identify loci used to bootstrap TCN means...done Number of (#hets, #homs, #nonSNPs): (2108,0,5478) Bootstrapping while preserving (#hets, #homs, #nonSNPs)... Number of bootstrap samples: 100 Bootstrapping while preserving (#hets, #homs, #nonSNPs)...done Segment #1 (chr 1, tcnId=1, dhId=1) of 4...done Segment #2 (chr 1, tcnId=2, dhId=1) of 4... chromosome tcnId dhId tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs 2 1 2 1 120992604 141510002 0 NA 0 tcnNbrOfHets dhStart dhEnd dhNbrOfLoci dhMean c1Mean c2Mean 2 0 NA NA 0 NA NA NA Number of TCNs: 0 Number of DHs: 0 int 0 int [1:14658] 1 2 3 4 5 6 7 8 9 10 ... int(0) Identify loci used to bootstrap DH means... Heterozygous SNPs to resample for DH: int 0 Identify loci used to bootstrap DH means...done Identify loci used to bootstrap TCN means... SNPs: int(0) Non-polymorphic loci: int(0) Heterozygous SNPs to resample for TCN: int(0) Homozygous SNPs to resample for TCN: int(0) Non-polymorphic loci to resample for TCN: int(0) Heterozygous SNPs with non-DH to resample for TCN: int(0) Loci to resample for TCN: int(0) Identify loci used to bootstrap TCN means...done Number of (#hets, #homs, #nonSNPs): (0,0,0) Bootstrapping while preserving (#hets, #homs, #nonSNPs)... Number of bootstrap samples: 100 Bootstrapping while preserving (#hets, #homs, #nonSNPs)...done Segment #2 (chr 1, tcnId=2, dhId=1) of 4...done Segment #3 (chr 1, tcnId=3, dhId=1) of 4... chromosome tcnId dhId tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs 3 1 3 1 141510003 185449813 2681 2.0689 777 tcnNbrOfHets dhStart dhEnd dhNbrOfLoci dhMean c1Mean c2Mean 3 777 141510003 185449813 777 0.0973 0.933798 1.135102 Number of TCNs: 2681 Number of DHs: 777 int [1:2681] 7587 7588 7589 7590 7591 7592 7593 7594 7595 7596 ... int [1:14658] 1 2 3 4 5 6 7 8 9 10 ... int [1:2681] 7587 7588 7589 7590 7591 7592 7593 7594 7595 7596 ... Identify loci used to bootstrap DH means... Heterozygous SNPs to resample for DH: int [1:777] 7587 7588 7594 7614 7616 7626 7627 7628 7635 7638 ... Identify loci used to bootstrap DH means...done Identify loci used to bootstrap TCN means... SNPs: int [1:777] 7587 7588 7594 7614 7616 7626 7627 7628 7635 7638 ... Non-polymorphic loci: int [1:1904] 7589 7590 7591 7592 7593 7595 7596 7597 7598 7599 ... Heterozygous SNPs to resample for TCN: int [1:777] 7587 7588 7594 7614 7616 7626 7627 7628 7635 7638 ... Homozygous SNPs to resample for TCN: int(0) Non-polymorphic loci to resample for TCN: int [1:1904] 7589 7590 7591 7592 7593 7595 7596 7597 7598 7599 ... Heterozygous SNPs with non-DH to resample for TCN: int(0) Loci to resample for TCN: int [1:2681] 7587 7588 7589 7590 7591 7592 7593 7594 7595 7596 ... Identify loci used to bootstrap TCN means...done Number of (#hets, #homs, #nonSNPs): (777,0,1904) Bootstrapping while preserving (#hets, #homs, #nonSNPs)... Number of bootstrap samples: 100 Bootstrapping while preserving (#hets, #homs, #nonSNPs)...done Segment #3 (chr 1, tcnId=3, dhId=1) of 4...done Segment #4 (chr 1, tcnId=4, dhId=1) of 4... chromosome tcnId dhId tcnStart tcnEnd tcnNbrOfLoci tcnMean tcnNbrOfSNPs 4 1 4 1 185449813 247137334 4391 2.6341 1311 tcnNbrOfHets dhStart dhEnd dhNbrOfLoci dhMean c1Mean c2Mean 4 1311 185449813 247137334 1311 0.2295 1.014787 1.619313 Number of TCNs: 4391 Number of DHs: 1311 int [1:4391] 10268 10269 10270 10271 10272 10273 10274 10275 10276 10277 ... int [1:14658] 1 2 3 4 5 6 7 8 9 10 ... int [1:4391] 10268 10269 10270 10271 10272 10273 10274 10275 10276 10277 ... Identify loci used to bootstrap DH means... Heterozygous SNPs to resample for DH: int [1:1311] 10269 10271 10281 10282 10283 10293 10295 10297 10300 10302 ... Identify loci used to bootstrap DH means...done Identify loci used to bootstrap TCN means... SNPs: int [1:1311] 10269 10271 10281 10282 10283 10293 10295 10297 10300 10302 ... Non-polymorphic loci: int [1:3080] 10268 10270 10272 10273 10274 10275 10276 10277 10278 10279 ... Heterozygous SNPs to resample for TCN: int [1:1311] 10269 10271 10281 10282 10283 10293 10295 10297 10300 10302 ... Homozygous SNPs to resample for TCN: int(0) Non-polymorphic loci to resample for TCN: int [1:3080] 10268 10270 10272 10273 10274 10275 10276 10277 10278 10279 ... Heterozygous SNPs with non-DH to resample for TCN: int(0) Loci to resample for TCN: int [1:4391] 10268 10269 10270 10271 10272 10273 10274 10275 10276 10277 ... Identify loci used to bootstrap TCN means...done Number of (#hets, #homs, #nonSNPs): (1311,0,3080) Bootstrapping while preserving (#hets, #homs, #nonSNPs)... Number of bootstrap samples: 100 Bootstrapping while preserving (#hets, #homs, #nonSNPs)...done Segment #4 (chr 1, tcnId=4, dhId=1) of 4...done Bootstrapped segment mean levels num [1:4, 1:100, 1:4] 1.39 NA 2.08 2.63 1.38 ... - attr(*, "dimnames")=List of 3 ..$ : NULL ..$ : NULL ..$ : chr [1:4] "tcn" "dh" "c1" "c2" Calculating (C1,C2) mean levels from (TCN,DH) mean levels... num [1:4, 1:100, 1:4] 1.39 NA 2.08 2.63 1.38 ... - attr(*, "dimnames")=List of 3 ..$ : NULL ..$ : NULL ..$ : chr [1:4] "tcn" "dh" "c1" "c2" Calculating (C1,C2) mean levels from (TCN,DH) mean levels...done Calculating polar (alpha,radius,manhattan) for change points... num [1:3, 1:100, 1:2] NA NA -0.0752 NA NA ... - attr(*, "dimnames")=List of 3 ..$ : NULL ..$ : NULL ..$ : chr [1:2] "c1" "c2" Bootstrapped change points num [1:3, 1:100, 1:5] NA NA -1.73 NA NA ... - attr(*, "dimnames")=List of 3 ..$ : NULL ..$ : NULL ..$ : chr [1:5] "alpha" "radius" "manhattan" "d1" ... Calculating polar (alpha,radius,manhattan) for change points...done Bootstrapping (TCN,DH,C1,C2) segment mean levels...done Summarizing bootstrapped segment ('tcn', 'dh', 'c1', 'c2') data... num [1:4, 1:4, 1:4] NA NA NA NA NA NA NA NA NA NA ... - attr(*, "dimnames")=List of 3 ..$ : NULL ..$ : chr [1:4] "2.5%" "5%" "95%" "97.5%" ..$ : chr [1:4] "tcn" "dh" "c1" "c2" Field #1 ('tcn') of 4... Segment #1 of 4... Segment #1 of 4...done Segment #2 of 4... Segment #2 of 4...done Segment #3 of 4... Segment #3 of 4...done Segment #4 of 4... Segment #4 of 4...done Field #1 ('tcn') of 4...done Field #2 ('dh') of 4... Segment #1 of 4... Segment #1 of 4...done Segment #2 of 4... Segment #2 of 4...done Segment #3 of 4... Segment #3 of 4...done Segment #4 of 4... Segment #4 of 4...done Field #2 ('dh') of 4...done Field #3 ('c1') of 4... Segment #1 of 4... Segment #1 of 4...done Segment #2 of 4... Segment #2 of 4...done Segment #3 of 4... Segment #3 of 4...done Segment #4 of 4... Segment #4 of 4...done Field #3 ('c1') of 4...done Field #4 ('c2') of 4... Segment #1 of 4... Segment #1 of 4...done Segment #2 of 4... Segment #2 of 4...done Segment #3 of 4... Segment #3 of 4...done Segment #4 of 4... Segment #4 of 4...done Field #4 ('c2') of 4...done Bootstrap statistics num [1:4, 1:4, 1:4] 1.38 NA 2.06 2.63 1.38 ... - attr(*, "dimnames")=List of 3 ..$ : NULL ..$ : chr [1:4] "2.5%" "5%" "95%" "97.5%" ..$ : chr [1:4] "tcn" "dh" "c1" "c2" Statistical sanity checks (iff B >= 100)... Available summaries: 2.5%, 5%, 95%, 97.5% Available quantiles: 0.025, 0.05, 0.95, 0.975 num [1:4, 1:4, 1:4] 1.38 NA 2.06 2.63 1.38 ... - attr(*, "dimnames")=List of 3 ..$ : NULL ..$ : chr [1:4] "2.5%" "5%" "95%" "97.5%" ..$ : chr [1:4] "tcn" "dh" "c1" "c2" Field #1 ('tcn') of 4... Seg 1. mean=1.3853, range=[1.37909,1.39287], n=7586 Seg 2. mean=NA, range=[NA,NA], n=0 Seg 3. mean=2.0689, range=[2.05903,2.079], n=2681 Seg 4. mean=2.6341, range=[2.62504,2.64649], n=4391 Field #1 ('tcn') of 4...done Field #2 ('dh') of 4... Seg 1. mean=0.5116, range=[0.502148,0.519941], n=2108 Seg 2. mean=NA, range=[NA,NA], n=NA Seg 3. mean=0.0973, range=[0.0906366,0.105818], n=777 Seg 4. mean=0.2295, range=[0.222919,0.237005], n=1311 Field #2 ('dh') of 4...done Field #3 ('c1') of 4... Seg 1. mean=0.33829, range=[0.332209,0.345936], n=2108 Seg 2. mean=NA, range=[NA,NA], n=NA Seg 3. mean=0.933798, range=[0.924112,0.941776], n=777 Seg 4. mean=1.01479, range=[1.00381,1.02461], n=1311 Field #3 ('c1') of 4...done Field #4 ('c2') of 4... Seg 1. mean=1.04701, range=[1.03882,1.05318], n=2108 Seg 2. mean=NA, range=[NA,NA], n=NA Seg 3. mean=1.1351, range=[1.12454,1.1465], n=777 Seg 4. mean=1.61931, range=[1.60862,1.63328], n=1311 Field #4 ('c2') of 4...done Statistical sanity checks (iff B >= 100)...done Summarizing bootstrapped segment ('tcn', 'dh', 'c1', 'c2') data...done Summarizing bootstrapped changepoint ('alpha', 'radius', 'manhattan', 'd1', 'd2') data... num [1:3, 1:4, 1:5] NA NA NA NA NA NA NA NA NA NA ... - attr(*, "dimnames")=List of 3 ..$ : NULL ..$ : chr [1:4] "2.5%" "5%" "95%" "97.5%" ..$ : chr [1:5] "alpha" "radius" "manhattan" "d1" ... Field #1 ('alpha') of 5... Changepoint #1 of 3... Changepoint #1 of 3...done Changepoint #2 of 3... Changepoint #2 of 3...done Changepoint #3 of 3... Changepoint #3 of 3...done Field #1 ('alpha') of 5...done Field #2 ('radius') of 5... Changepoint #1 of 3... Changepoint #1 of 3...done Changepoint #2 of 3... Changepoint #2 of 3...done Changepoint #3 of 3... Changepoint #3 of 3...done Field #2 ('radius') of 5...done Field #3 ('manhattan') of 5... Changepoint #1 of 3... Changepoint #1 of 3...done Changepoint #2 of 3... Changepoint #2 of 3...done Changepoint #3 of 3... Changepoint #3 of 3...done Field #3 ('manhattan') of 5...done Field #4 ('d1') of 5... Changepoint #1 of 3... Changepoint #1 of 3...done Changepoint #2 of 3... Changepoint #2 of 3...done Changepoint #3 of 3... Changepoint #3 of 3...done Field #4 ('d1') of 5...done Field #5 ('d2') of 5... Changepoint #1 of 3... Changepoint #1 of 3...done Changepoint #2 of 3... Changepoint #2 of 3...done Changepoint #3 of 3... Changepoint #3 of 3...done Field #5 ('d2') of 5...done Bootstrap statistics num [1:3, 1:4, 1:5] NA NA -1.77 NA NA ... - attr(*, "dimnames")=List of 3 ..$ : NULL ..$ : chr [1:4] "2.5%" "5%" "95%" "97.5%" ..$ : chr [1:5] "alpha" "radius" "manhattan" "d1" ... Summarizing bootstrapped changepoint ('alpha', 'radius', 'manhattan', 'd1', 'd2') data...done Resample (TCN,DH) signals and re-estimate summaries for segment & changepoint...done > print(fit) chromosome tcnId dhId start end tcnNbrOfLoci tcnMean tcnNbrOfSNPs 1 1 1 1 554484 120992603 7586 1.3853 2108 2 1 2 1 120992604 141510002 0 NA 0 3 1 3 1 141510003 185449813 2681 2.0689 777 4 1 4 1 185449813 247137334 4391 2.6341 1311 tcnNbrOfHets dhNbrOfLoci dhMean c1Mean c2Mean 1 2108 2108 0.5116 0.3382903 1.047010 2 0 NA NA NA NA 3 777 777 0.0973 0.9337980 1.135102 4 1311 1311 0.2295 1.0147870 1.619313 > plotTracks(fit) > > > > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > # Calling segments with run of homozygosity (ROH) > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > fit <- callROH(fit, verbose=-10) Calling ROH... Segment #1 of 4... Calling ROH for a single segment... Number of SNPs: 7586 Calling ROH for a single segment...done Segment #1 of 4...done Segment #2 of 4... Calling ROH for a single segment... Number of SNPs: 0 Calling ROH for a single segment...done Segment #2 of 4...done Segment #3 of 4... Calling ROH for a single segment... Number of SNPs: 2681 Calling ROH for a single segment...done Segment #3 of 4...done Segment #4 of 4... Calling ROH for a single segment... Number of SNPs: 4391 Calling ROH for a single segment...done Segment #4 of 4...done ROH calls: logi [1:4] FALSE NA FALSE FALSE Mode FALSE NA's logical 3 1 Calling ROH...done > print(fit) chromosome tcnId dhId start end tcnNbrOfLoci tcnMean tcnNbrOfSNPs 1 1 1 1 554484 120992603 7586 1.3853 2108 2 1 2 1 120992604 141510002 0 NA 0 3 1 3 1 141510003 185449813 2681 2.0689 777 4 1 4 1 185449813 247137334 4391 2.6341 1311 tcnNbrOfHets dhNbrOfLoci dhMean c1Mean c2Mean rohCall 1 2108 2108 0.5116 0.3382903 1.047010 FALSE 2 0 NA NA NA NA NA 3 777 777 0.0973 0.9337980 1.135102 FALSE 4 1311 1311 0.2295 1.0147870 1.619313 FALSE > plotTracks(fit) > > > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > # Estimate background > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > kappa <- estimateKappa(fit, verbose=-10) Estimate global background (including normal contamination and more)... Number of segments: 3 Estimating threshold Delta0.5 from the empirical density of C1:s... adjust: 1 minDensity: 0.2 ploidy: 2 All peaks: type x density 1 peak 0.3362194 1.101242 3 peak 0.9811492 1.065635 C1=0 and C1=1 peaks: type x density 1 peak 0.3362194 1.101242 3 peak 0.9811492 1.065635 Estimate of Delta0.5: 0.65868427808456 Estimating threshold Delta0.5 from the empirical density of C1:s...done Number of segments with C1 < Delta0.5: 1 Estimate of kappa: 0.33829026 Estimate global background (including normal contamination and more)...done Warning message: In density.default(c1, weights = weights, adjust = adjust, from = from, : Selecting bandwidth *not* using 'weights' > print(kappa) [1] 0.3382903 > ## [1] 0.226011 > > > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > # Calling segments in allelic balance (AB) > # NOTE: Ideally, this should be done on whole-genome data > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > # Explicitly estimate the threshold in DH for calling AB > # (which be done by default by the caller, if skipped here) > deltaAB <- estimateDeltaAB(fit, flavor="qq(DH)", verbose=-10) Estimating DH threshold for calling allelic imbalances... flavor: qq(DH) scale: 1 Estimating DH threshold for AB caller... quantile #1: 0.05 Symmetric quantile #2: 0.9 Number of segments: 3 Weighted 5% quantile of DH: 0.257710 Number of segments with small DH: 2 Number of data points: 7072 Number of finite data points: 2088 Estimate of (1-0.9):th and 50% quantiles: (0.0310411,0.163658) Estimate of 0.9:th "symmetric" quantile: 0.296275 Estimating DH threshold for AB caller...done Estimated delta: 0.296 Estimating DH threshold for calling allelic imbalances...done > if (Sys.getenv("_R_CHECK_FULL_") == "") { + # Ad hoc workaround for not utilizing all of the data + # in the test, which results in a poor estimate + deltaAB <- 0.165 + } > print(deltaAB) [1] 0.165 > ## [1] 0.1657131 > > fit <- callAB(fit, delta=deltaAB, verbose=-10) Calling segments of allelic balance from one-sided DH bootstrap confidence intervals... delta (offset adjusting for bias in DH): 0.165 alpha (CI quantile; significance level): 0.05 Calling segments... Number of segments called allelic balance (AB): 1 (25.00%) of 4 Calling segments...done Calling segments of allelic balance from one-sided DH bootstrap confidence intervals...done > print(fit) chromosome tcnId dhId start end tcnNbrOfLoci tcnMean tcnNbrOfSNPs 1 1 1 1 554484 120992603 7586 1.3853 2108 2 1 2 1 120992604 141510002 0 NA 0 3 1 3 1 141510003 185449813 2681 2.0689 777 4 1 4 1 185449813 247137334 4391 2.6341 1311 tcnNbrOfHets dhNbrOfLoci dhMean c1Mean c2Mean rohCall abCall 1 2108 2108 0.5116 0.3382903 1.047010 FALSE FALSE 2 0 NA NA NA NA NA NA 3 777 777 0.0973 0.9337980 1.135102 FALSE TRUE 4 1311 1311 0.2295 1.0147870 1.619313 FALSE FALSE > plotTracks(fit) > > # Even if not explicitly specified, the estimated > # threshold parameter is returned by the caller > stopifnot(fit$params$deltaAB == deltaAB) > > > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > # Calling segments in loss-of-heterozygosity (LOH) > # NOTE: Ideally, this should be done on whole-genome data > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > # Explicitly estimate the threshold in C1 for calling LOH > # (which be done by default by the caller, if skipped here) > deltaLOH <- estimateDeltaLOH(fit, flavor="minC1|nonAB", verbose=-10) Estimating DH threshold for calling LOH... flavor: minC1|nonAB Estimating DH threshold for calling LOH as the midpoint between guessed C1=0 and C1=1... Argument 'midpoint': 0.5 Number of segments: 4 Number of segments in allelic balance: 1 (25.0%) of 4 Number of segments not in allelic balance: 2 (50.0%) of 4 Number of segments in allelic balance and TCN <= 3.00: 1 (25.0%) of 4 C: 2.07 Corrected C1 (=C/2): 1.03 Number of DHs: 777 Weights: 1 Weighted median of (corrected) C1 in allelic balance: 1.034 Smallest C1 among segments not in allelic balance: 0.338 There are 1 segments with in total 2108 heterozygous SNPs with this level. Midpoint between the two: 0.686 Estimating DH threshold for calling LOH as the midpoint between guessed C1=0 and C1=1...done delta: 0.686 Estimating DH threshold for calling LOH...done > print(deltaLOH) [1] 0.6863701 > ## [1] 0.625175 > > fit <- callLOH(fit, delta=deltaLOH, verbose=-10) Calling segments of allelic balance from one-sided DH bootstrap confidence intervals... delta (offset adjusting for bias in C1): 0.68637013 alpha (CI quantile; significance level): 0.05 Calling segments... Number of segments called low C1 (LowC1, "LOH_C1"): 1 (25.00%) of 4 Calling segments...done Calling segments of allelic balance from one-sided DH bootstrap confidence intervals...done > print(fit) chromosome tcnId dhId start end tcnNbrOfLoci tcnMean tcnNbrOfSNPs 1 1 1 1 554484 120992603 7586 1.3853 2108 2 1 2 1 120992604 141510002 0 NA 0 3 1 3 1 141510003 185449813 2681 2.0689 777 4 1 4 1 185449813 247137334 4391 2.6341 1311 tcnNbrOfHets dhNbrOfLoci dhMean c1Mean c2Mean rohCall abCall lohCall 1 2108 2108 0.5116 0.3382903 1.047010 FALSE FALSE TRUE 2 0 NA NA NA NA NA NA NA 3 777 777 0.0973 0.9337980 1.135102 FALSE TRUE FALSE 4 1311 1311 0.2295 1.0147870 1.619313 FALSE FALSE FALSE > plotTracks(fit) > > # Even if not explicitly specified, the estimated > # threshold parameter is returned by the caller > stopifnot(fit$params$deltaLOH == deltaLOH) > > > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > # Calling segments that are gained, copy neutral, and lost > # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - > fit <- callGNL(fit, verbose=-10) Calling gain, neutral, and loss based TCNs of AB segments... Calling neutral TCNs... callCopyNeutralByTCNofAB... Alpha: 0.05 Delta CN: 0.33085487 Calling copy-neutral segments... Retrieve TCN confidence intervals for all segments... Interval: [0.025,0.975] Retrieve TCN confidence intervals for all segments...done Estimating TCN confidence interval of copy-neutral AB segments... calcStatsForCopyNeutralABs... Identifying copy neutral AB segments... Number of AB segments: 1 Identifying segments that are copy neutral states... Number of segments in allelic balance: 1 Identifying segments that are copy neutral states...done Number of copy-neutral AB segments: 1 Extracting all copy neutral AB segments across all chromosomes into one big segment... chromosome tcnId dhId start end tcnNbrOfLoci tcnMean tcnNbrOfSNPs 3 1 3 1 141510003 185449813 2681 2.0689 777 tcnNbrOfHets dhNbrOfLoci dhMean c1Mean c2Mean rohCall abCall lohCall 3 777 777 0.0973 0.933798 1.135102 FALSE TRUE FALSE Extracting all copy neutral AB segments across all chromosomes into one big segment...done Identifying copy neutral AB segments...done Bootstrap the identified copy-neutral states... Bootstrap the identified copy-neutral states...done tcnNbrOfLoci tcnMean tcnNbrOfSNPs tcnNbrOfHets dhNbrOfLoci dhMean c1Mean 3 2681 2.0689 777 777 777 0.0973 0.933798 c2Mean tcn_2.5% tcn_5% tcn_95% tcn_97.5% dh_2.5% dh_5% dh_95% 3 1.135102 2.055164 2.057694 2.078831 2.081454 0.08974138 0.09080508 0.1035891 dh_97.5% c1_2.5% c1_5% c1_95% c1_97.5% c2_2.5% c2_5% c2_95% 3 0.1050478 0.923788 0.925412 0.9417056 0.9433752 1.124908 1.126631 1.143571 c2_97.5% 3 1.145214 calcStatsForCopyNeutralABs...done Bootstrap statistics for copy-neutral AB segments: tcnNbrOfLoci tcnMean tcnNbrOfSNPs tcnNbrOfHets dhNbrOfLoci dhMean c1Mean 3 2681 2.0689 777 777 777 0.0973 0.933798 c2Mean tcn_2.5% tcn_5% tcn_95% tcn_97.5% dh_2.5% dh_5% dh_95% 3 1.135102 2.055164 2.057694 2.078831 2.081454 0.08974138 0.09080508 0.1035891 dh_97.5% c1_2.5% c1_5% c1_95% c1_97.5% c2_2.5% c2_5% c2_95% 3 0.1050478 0.923788 0.925412 0.9417056 0.9433752 1.124908 1.126631 1.143571 c2_97.5% 3 1.145214 [1] "TCN statistics:" tcnMean tcn_2.5% tcn_5% tcn_95% tcn_97.5% 2.068900 2.055164 2.057694 2.078831 2.081454 95%-confidence interval of TCN mean for the copy-neutral state: [2.05516,2.08145] (mean=2.0689) Estimating TCN confidence interval of copy-neutral AB segments...done Identify all copy-neutral segments... DeltaCN: +/-0.330855 Call ("acceptance") region: [1.72431,2.41231] Total number of segments: 4 Number of segments called allelic balance: 1 Number of segments called copy neutral: 1 Number of AB segments called copy neutral: 1 Number of non-AB segments called copy neutral: 0 Identify all copy-neutral segments...done Calling copy-neutral segments...done callCopyNeutralByTCNofAB...done Calling neutral TCNs...done Number of NTCN calls: 1 (25.00%) of 4 Mean TCN of AB segments: 2.06831 Calling loss... Number of loss calls: 1 (25.00%) of 4 Calling loss...done Calling gain... Number of loss calls: 1 (25.00%) of 4 Calling gain...done Calling gain, neutral, and loss based TCNs of AB segments...done Warning message: In density.default(c1, weights = weights, adjust = adjust, from = from, : Selecting bandwidth *not* using 'weights' > print(fit) chromosome tcnId dhId start end tcnNbrOfLoci tcnMean tcnNbrOfSNPs 1 1 1 1 554484 120992603 7586 1.3853 2108 2 1 2 1 120992604 141510002 0 NA 0 3 1 3 1 141510003 185449813 2681 2.0689 777 4 1 4 1 185449813 247137334 4391 2.6341 1311 tcnNbrOfHets dhNbrOfLoci dhMean c1Mean c2Mean rohCall abCall lohCall 1 2108 2108 0.5116 0.3382903 1.047010 FALSE FALSE TRUE 2 0 NA NA NA NA NA NA NA 3 777 777 0.0973 0.9337980 1.135102 FALSE TRUE FALSE 4 1311 1311 0.2295 1.0147870 1.619313 FALSE FALSE FALSE ntcnCall lossCall gainCall 1 FALSE TRUE FALSE 2 NA NA NA 3 TRUE FALSE FALSE 4 FALSE FALSE TRUE > plotTracks(fit) > > proc.time() user system elapsed 4.03 0.20 4.79