seeking help on cghMCR pleeeease
0
0
Entering edit mode
nac ▴ 280
@nac-4545
Last seen 10.3 years ago
HI, I am using cghMCR and I would need some advice on the software behaviour. I don't seem to get any different results when I change the gapAllowed parameter. I have used DNAcopy to create a segmented data for all my samples. Then I choose to use cghMCR to get common alterations between my samples. These are all arguments from the package doc(latest version): segments is a data frame extracted from the "output" element of the object returned by segment of the package DNAcopy or getSegments gapAllowed is an integer specifying low threshold of base pair number to separate two adjacent segments, belower which the two segments will be joined as an altered span alteredLow is a positive number between 0 and 1 specifying the lower resh- old percential value. Only segments with values falling below this threshold are considered as altered span alteredHigh is a positive number between 0 and 1 specifying the upper resh- old percential value. Only segments with values falling over this threshold are considered as altered span recurrence is an integer between 1 and 100 that specifies the rate of occur- rence for a gain or loss that are observed across sample. Only gains/losses with ocurrence rate grater than the threshold values are declared as MCRs spanLimit is an integer that defines the leangh of altered spans that can be considered as locus. It is not of any use at this time thresholdType is a character string that can be either "quantile" or "value" indicating wether alteredLow or alteredHigh is quantial or actual value In my analysis, I have used cghMCR with threshold value log2R of low -0.25 and high =0.25 , a spanLimit of 2.10^7 and I have tried several gapAllowed values, 5, 500 , 5000 and 50000 .Each time I used the MCR function to identify the minimum regions. I have put an example of each results only on 25 lines of chromosome 4 to avoid massive sized files (see attached) , but basically whatever the gapAllowed size i get the same MCRs at the end for all postitions. I would expect this to vary as segments should be fused together differently given this parameter.. Could you please help with this and advise on what I might be doing wrong? Is spanLimit any use? in the doc, It is written "It is not of any use at this time"?? Another point, the gapAllowed is specified " is an integer specifying low threshold of base pair number " is the unit in base pair number, in several examples it seems that the units are in kb??? thanks a lot , code for gapAllowed=5 ##get the cghMCR function with these parameters cghmcr0.25T_5k_2_20M_sdundo1.5 <- cghMCR(sdundo.segData_order1.5, gapAllowed = 5, alteredLow = -0.25,alteredHigh = 0.25, spanLimit=20000000,recurrence=2,thresholdType=c("value")) ##identify the MCRs mcrs0.25T_5k_2_20M_sdundo1.5<-MCR(cghmcr0.25T_5k_2_20M_sdundo1.5) mcrs0.25T_5k_2_20M_sdundo.bind1.5 <-cbind ( mcrs0.25T_5k_2_20M_sdundo1.5[, c ( "chromosome", "status", "loc.start", "loc.end","mcr.start", "mcr.end","samples" ) ],z[1:9171] ) write.table(mcrs0.25T_5k_2_20M_sdundo.bind1.5, file="PT_mcrs0.25T_5k_2_20M_sdundo.bind1.5.txt", sep="\t", row.names=F) PT_mcrs0.25T_5k_2_20M_sdundo.bind1.5=read.table(file="PT_mcrs0.25T_5k_ 2_20M_sdundo.bind1.5.txt", sep="\t", header=T) ##only 25 lines of chromsome 4 test.5k=head(PT_mcrs0.25T_5k_2_20M_sdundo.bind1.5[PT_mcrs0.25T_5k_2_20 M_sdundo.bind1.5==4,],25) ##attached data write.table(test.5k, file="test.5k.txt", sep="\t") code for gapAllowed=500 cghmcr0.25T_500k_2_20M_sdundo1.5 <- cghMCR(sdundo.segData_order1.5, gapAllowed = 500, alteredLow = -0.25,alteredHigh = 0.25, spanLimit=20000000,recurrence=2,thresholdType=c("value")) mcrs0.25T_500k_2_20M_sdundo1.5<-MCR(cghmcr0.25T_500k_2_20M_sdundo1.5) mcrs0.25T_500k_2_20M_sdundo.bind1.5 <-cbind ( mcrs0.25T_500k_2_20M_sdundo1.5[, c ( "chromosome", "status", "loc.start", "loc.end","mcr.start", "mcr.end","samples" ) ],z[1:9171] ) write.table(mcrs0.25T_500k_2_20M_sdundo.bind1.5, file="PT_mcrs0.25T_500k_2_20M_sdundo.bind1.5.txt", sep="\t", row.names=F) PT_mcrs0.25T_500k_2_20M_sdundo.bind1.5=read.table(file="PT_mcrs0.25T_5 00k_2_20M_sdundo.bind1.5.txt", sep="\t", header=T) test.500k=head(PT_mcrs0.25T_500k_2_20M_sdundo.bind1.5[PT_mcrs0.25T_500 k_2_20M_sdundo.bind1.5==4,],25) write.table(test.500k, file="test.500k.txt", sep="\t") code for gapAllowed=5000 cghmcr0.25T_5000k_2_20M_sdundo1.5 <- cghMCR(sdundo.segData_order1.5, gapAllowed = 5000, alteredLow = -0.25,alteredHigh = 0.25, spanLimit=20000000,recurrence=2,thresholdType=c("value")) mcrs0.25T_5000k_2_20M_sdundo1.5<-MCR(cghmcr0.25T_5000k_2_20M_sdundo1.5 ) mcrs0.25T_5000k_2_20M_sdundo.bind1.5 <-cbind ( mcrs0.25T_5000k_2_20M_sdundo1.5[, c ( "chromosome", "status", "loc.start", "loc.end","mcr.start", "mcr.end","samples" ) ],z[1:9171] ) write.table(mcrs0.25T_5000k_2_20M_sdundo.bind1.5, file="PT_mcrs0.25T_5000k_2_20M_sdundo.bind1.5.txt", sep="\t", row.names=F) PT_mcrs0.25T_5000k_2_20M_sdundo.bind1.5=read.table(file="PT_mcrs0.25T_ 5000k_2_20M_sdundo.bind1.5.txt", sep="\t", header=T) test.5000k=head(PT_mcrs0.25T_5000k_2_20M_sdundo.bind1.5[PT_mcrs0.25T_5 000k_2_20M_sdundo.bind1.5==4,],25) write.table(test.5000k, file="test.5000k.txt", sep="\t") code for gapAllowed=500000 cghmcr0.25T_500000k_2_20M_sdundo1.5 <- cghMCR(sdundo.segData_order1.5, gapAllowed = 500000, alteredLow = -0.25,alteredHigh = 0.25, spanLimit=20000000,recurrence=2,thresholdType=c("value")) mcrs0.25T_500000k_2_20M_sdundo1.5<-MCR(cghmcr0.25T_500000k_2_20M_sdund o1.5) mcrs0.25T_500000k_2_20M_sdundo.bind1.5 <-cbind ( mcrs0.25T_500000k_2_20M_sdundo1.5[, c ( "chromosome", "status", "loc.start", "loc.end","mcr.start", "mcr.end","samples" ) ],z[1:9171] ) write.table(mcrs0.25T_500000k_2_20M_sdundo.bind1.5, file="PT_mcrs0.25T_500000k_2_20M_sdundo.bind1.5.txt", sep="\t", row.names=F) PT_mcrs0.25T_500000k_2_20M_sdundo.bind1.5=read.table(file="PT_mcrs0.25 T_500000k_2_20M_sdundo.bind1.5.txt", sep="\t", header=T) test.500000k=head(PT_mcrs0.25T_500000k_2_20M_sdundo.bind1.5[PT_mcrs0.2 5T_500000k_2_20M_sdundo.bind1.5==4,],25) write.table(test.500000k, file="test.500000k.txt", sep="\t") sessioninfo() R version 2.13.0 (2011-04-13) Platform: x86_64-unknown-linux-gnu (64-bit) locale: [1] LC_CTYPE=en_GB.UTF-8 LC_NUMERIC=C [3] LC_TIME=en_GB.UTF-8 LC_COLLATE=en_GB.UTF-8 [5] LC_MONETARY=C LC_MESSAGES=C [7] LC_PAPER=en_GB.UTF-8 LC_NAME=C [9] LC_ADDRESS=C LC_TELEPHONE=C [11] LC_MEASUREMENT=en_GB.UTF-8 LC_IDENTIFICATION=C attached base packages: [1] tools stats graphics grDevices utils datasets methods [8] base other attached packages: [1] cghMCR_1.10.0 limma_3.4.5 CNTools_1.6.0 genefilter_1.30.0 [5] DNAcopy_1.22.1 loaded via a namespace (and not attached): [1] annotate_1.26.1 AnnotationDbi_1.10.1 Biobase_2.8.0 [4] DBI_0.2-5 RSQLite_0.9-1 splines_2.13.0 [7] survival_2.35-8 xtable_1.6-0 -- The Wellcome Trust Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is 215 Euston Road, London, NW1 2BE. -------------- next part -------------- An embedded and charset-unspecified text was scrubbed... Name: test.5k.txt URL: <https: stat.ethz.ch="" pipermail="" bioconductor="" attachments="" 20111212="" 2c70ceda="" attachment.txt=""> -------------- next part -------------- An embedded and charset-unspecified text was scrubbed... Name: test.500k.txt URL: <https: stat.ethz.ch="" pipermail="" bioconductor="" attachments="" 20111212="" 2c70ceda="" attachment-0001.txt=""> -------------- next part -------------- An embedded and charset-unspecified text was scrubbed... Name: test.5000k.txt URL: <https: stat.ethz.ch="" pipermail="" bioconductor="" attachments="" 20111212="" 2c70ceda="" attachment-0002.txt=""> -------------- next part -------------- An embedded and charset-unspecified text was scrubbed... Name: test.500000k.txt URL: <https: stat.ethz.ch="" pipermail="" bioconductor="" attachments="" 20111212="" 2c70ceda="" attachment-0003.txt="">
DNAcopy cghMCR DNAcopy cghMCR • 1.1k views
ADD COMMENT

Login before adding your answer.

Traffic: 517 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6