Hi, with help of the data you sent me, I could reproduce the issue, identify and correct the bug. Update to PSCBS v0.19.4:
http://cran.r-project.org/web/packages/PSCBS/ Let me know if it solved it for you /Henrik On Fri, Jan 6, 2012 at 5:43 AM, Henrik Bengtsson <h...@biostat.ucsf.edu> wrote: > Hi Minya, > > what does sessionInfo() report after library("PSCBS")? > > Also, could you please send me (offline; not to the mailinglist) the > following 'MP20120106.Rbin' file: > >> df <- data.frame(chromosome=chromosome, x=x, CT=CT, betaT=betaT, >> betaN=betaN); >> saveObject(df, "MP20120106.Rbin"); > > /Henrik > > > On Friday, January 6, 2012, Minya Pu wrote: >> >> Hi, >> >> I am doing paired PSCBS analysis on a Genome Wide SNP array 6.0 data >> set. Following the code on the aroma web site, among 11 pairs of the samples >> from, I have successfully done the analysis for 10 of them but got stuck in >> one. >> >> Below is the error message from the ‘bad’ pair, it came out at the >> segmentByPairedPSCBS step: >> >> >> >> startRow endRow >> >> 1 1 3876 >> >> 2 3877 8743 >> >> 3 8744 11707 >> >> 4 11708 11901 >> >> 5 11902 13071 >> >> 6 13072 13474 >> >> 7 13475 14880 >> >> 8 14881 16114 >> >> 9 16115 21428 >> >> 10 21429 21466 >> >> 11 21467 23886 >> >> 12 23887 24009 >> >> 13 24010 35079 >> >> startRow endRow >> >> 1 25 3853 >> >> 2 3896 8716 >> >> 3 8762 11707 >> >> 4 11710 11897 >> >> 5 11909 13068 >> >> 6 13075 13441 >> >> 7 13492 14874 >> >> 8 14883 16113 >> >> 9 16116 21428 >> >> 10 21430 21465 >> >> 11 21469 23879 >> >> 12 23890 23988 >> >> 13 24028 35025 >> >> 20120105 16:43:25| Number of TCNs before: 35078 >> >> 20120105 16:43:25| Number of TCNs after: 35079 >> >> 20120105 16:43:25| TCN segment #1 ('1') of 5...done >> >> 20120105 16:43:25| TCN segment #2 ('2') of 5... >> >> 20120105 16:43:25| Nothing todo. Only one DH segmentation. Skipping. >> >> 20120105 16:43:25| TCN segment #2 ('2') of 5...done >> >> 20120105 16:43:25| TCN segment #3 ('3') of 5... >> >> 20120105 16:43:25| Nothing todo. Only one DH segmentation. Skipping. >> >> 20120105 16:43:25| TCN segment #3 ('3') of 5...done >> >> 20120105 16:43:25| TCN segment #4 ('4') of 5... >> >> 20120105 16:43:25| Nothing todo. Only one DH segmentation. Skipping. >> >> 20120105 16:43:25| TCN segment #4 ('4') of 5...done >> >> 20120105 16:43:25| TCN segment #5 ('5') of 5... >> >> 20120105 16:43:25| Nothing todo. Only one DH segmentation. Skipping. >> >> 20120105 16:43:25| TCN segment #5 ('5') of 5...done >> >> Error: all(tcnSegRowsCC[-nrow(tcnSegRowsCC), 2] < tcnSegRowsCC[-1, 1], >> .... is not TRUE >> >> >> >> Any suggestion on what might go wrong? >> >> >> >> Below is the code I used. >> >> Thanks so much, >> >> Minya >> >> >> >> >> >> library("PSCBS"); >> >> >> >> for (i in 1:11 ){ >> >> >> >> # Extract tumor-normal pair of interest >> >> pair <- c(T=tumors[i], N=normals[i]); >> >> csR2 <- extract(csR, indexOf(csR, pair)); >> >> print(csR2); >> >> res <- doASCRMAv2(csR2, verbose=verbose); >> >> >> >> data <- extractPSCNArray(res$total); >> >> dimnames(data)[[3]] <- names(pair); >> >> str(data); >> >> # Total CNs for the tumor relative to the matched normal >> >> >> >> CT <- 2 * (data[,"total","T"] / data[,"total","N"]); >> >> >> >> # Allele B fractions for the tumor >> >> betaT <- data[,"fracB","T"]; >> >> # Allele B fractions for the normal >> >> betaN <- data[,"fracB","N"]; >> >> >> >> ugp <- getAromaUgpFile(res$total); >> >> >> >> chromosome <- ugp[,1,drop=TRUE]; >> >> >> >> x <- ugp[,2,drop=TRUE]; >> >> ### x is position ### >> >> >> >> df <- data.frame(chromosome=chromosome, x=x, CT=CT, betaT=betaT, >> betaN=betaN); >> >> >> >> fit1 <- segmentByPairedPSCBS(df, verbose=verbose); >> >> >> >> segs <- getSegments(fit1); >> >> print(segs); >> >> >> >> pairName =paste(pair, collapse="vs"); >> >> chrTag <- sprintf("Chr%s", seqToHumanReadable(getChromosomes(fit1))); >> >> toPNG(pairName, tags=c(chrTag, "PairedPSCBS"), width=840, aspectRatio=0.6, >> { >> >> plotTracks(fit1); >> >> }); >> >> >> >> } >> >> -- >> When reporting problems on aroma.affymetrix, make sure 1) to run the >> latest version of the package, 2) to report the output of sessionInfo() and >> traceback(), and 3) to post a complete code example. >> >> >> You received this message because you are subscribed to the Google Groups >> "aroma.affymetrix" group with website http://www.aroma-project.org/. >> To post to this group, send email to aroma-affymetrix@googlegroups.com >> To unsubscribe and other options, go to >> http://www.aroma-project.org/forum/ -- When reporting problems on aroma.affymetrix, make sure 1) to run the latest version of the package, 2) to report the output of sessionInfo() and traceback(), and 3) to post a complete code example. You received this message because you are subscribed to the Google Groups "aroma.affymetrix" group with website http://www.aroma-project.org/. To post to this group, send email to aroma-affymetrix@googlegroups.com To unsubscribe and other options, go to http://www.aroma-project.org/forum/