Hi Oscar,
I'm rather swamped right now, but just a few quick comments. Calling
so called copy-neutral segments via callNTCN() is a rather and
somewhat experimental method (I should state in the vignette). In
other words, it has not undergone equally much rigorous validation as
the other callers. Having said this, this caller is very of interest
to us (and also right now through several projects), so any feedback
is very useful. I will try to have a look at your data that you've
sent me offline, but I cannot promise when. From the error message,
it could simply be a bug rather than a problem with the
estimators/data.
About the avgDH="median". Exactly what commands did you call? Also,
you're using PSCBS v0.32.5, correct?
/Henrik
On Tue, Feb 12, 2013 at 8:06 AM, Oscar Rueda <oscar.ru...@cruk.cam.ac.uk> wrote:
> Hi Henrik,
>
> That works, but with another sample, trying the last step suggested in the
> vignette I get another error:
>
>> fit <- callLOH(fit, verbose = -10)
>> deltaCN <- estimateDeltaCN(fit, scale=1)
>> fit <- callNTCN(fit, delta=deltaCN, verbose=-10)
>
>
> ..........
> Identifying segments that are copy neutral states...done
> Number of copy-neutral AB segments: 2
> Extracting all copy neutral AB segments across all chromosomes into
> one big segment...
> chromosome tcnId dhId start end tcnNbrOfLoci tcnMean
> tcnNbrOfSNPs
> 128 4 21 1 70750718 70913014 83 1.4003
> 38
> 488 18 9 1 4105338 7613644 2635 1.5945
> 1546
> tcnNbrOfHets dhNbrOfLoci dhMean c1Mean c2Mean rohCall
> abCall lohCall
> 128 16 16 0.1261 0.6118611 0.7884389 FALSE
> TRUE NA
> 488 401 401 0.1307 0.6930494 0.9014506 FALSE
> TRUE NA
> Extracting all copy neutral AB segments across all chromosomes into
> one big segment...done
> Error: nbrOfSegments(fitNTCN) == 1 is not TRUE
> Identifying copy neutral AB segments...done
> calcStatsForCopyNeutralABs...done
> Estimating TCN confidence interval of copy-neutral AB segments...done
> Calling copy-neutral segments...done
> callCopyNeutralByTCNofAB...done
>
> By the way, I also tried
> avgDH="median"
> but I got this error
>
>
> Error in bootstrapTCNandDHByRegion.PairedPSCBS(fit, statsFcn = statsFcn, :
> INTERNAL ERROR: Incorrect DH mean!
>
>
>
>
> Thanks,
> Oscar
>
> Oscar M. Rueda, PhD.
> Postdoctoral Research Fellow, Caldas Lab, Breast Cancer Functional
> Genomics.
> University of Cambridge. Cancer Research UK Cambridge Institute.
> Li Ka Shing Centre, Robinson Way.
> Cambridge CB2 0RE
> England
>
> Please note that my email address will be changing to:
> oscar.ru...@cruk.cam.ac.uk
>
>
>
>>
>>
>>
>>
>>
>>
>>Henrik Bengtsson <henrik.bengts...@aroma-project.org>
>> Feb 09 01:45PM -0800
>>
>>
>>
>>
>>Hi Oscar,
>>
>>
>>
>>
>>
>>so you spotted a bug in the bootstrapping of the TCN mean-level
>>
>>
>>estimates. It only affected segments that had been called to be in
>>
>>
>>run-of-homozygosity (ROH) resulting in sampling from not all available
>>
>>
>>TCN loci. The effect was incorrect TCN mean quantile estimates for
>>
>>
>>those segments. For reasonable large segments this didn't matter
>>
>>
>>much, but for tiny segments (say < 5 loci) it made a small difference.
>>
>>
>>In your case the internal sanity checks caught this error in an ROH
>>
>>
>>segment of two loci.
>>
>>
>>
>>
>>
>>The bug has been fixed in PSCBS v0.32.5. You can grab it as usual:
>>
>>
>>
>>
>>
>>source("
>>http://aroma-project.org/hbLite.R
>>");
>>
>>
>>hbInstall("PSCBS");
>>
>>
>>
>>
>>
>>I've verified that it now works with the sample data you've sent me.
>>
>>
>>Let me know if it solves it for all 126-60=66 samples where you had
>>
>>
>>problems.
>>
>>
>>
>>
>>
>>
>>
>>
>>BTW, It is experimental/still early days, but you may want to move to
>>
>>
>>calculate DH using the median estimator instead of the default sample
>>
>>
>>mean. For more details, see 'Less biased Decrease of Heterozygosity
>>
>>
>>(DH) estimates' in Section 'Experimental' of the 'Paired PSCBS'
>>
>>
>>vignette of PSCBS v0.32.5. Also, if you have saved 'fit' objects, you
>>
>>
>>don't have to rerun the segmentation (=finding the change points) to
>>
>>
>>achieve this; you can do fit <- updateMeans(fit, avgDH="median") and
>>
>>
>>then redo the calling.
>>
>>
>>
>>
>>
>>
>>
>>
>>Thanks for reporting.
>>
>>
>>
>>
>>
>>Henrik
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>You received this message because you are subscribed to the Google Group
>>aroma-affymetrix.
>>You can
>>post via email <mailto:aroma-affymetrix@googlegroups.com>.
>>To unsubscribe from this group,
>>send <mailto:aroma-affymetrix+unsubscr...@googlegroups.com> an empty
>>message.
>>For more options,
>>visit <http://groups.google.com/group/aroma-affymetrix/topics> this group.
>>
>>--
>>--
>>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/ <http://www.aroma-project.org/forum/>
>>
>>---
>>You received this message because you are subscribed to the Google Groups
>>"aroma.affymetrix" group.
>>To unsubscribe from this group and stop receiving emails from it, send an
>>email to aroma-affymetrix+unsubscr...@googlegroups.com.
>>For more options, visit https://groups.google.com/groups/opt_out.
>>
>>
>>
>
> --
> --
> 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/
>
> ---
> You received this message because you are subscribed to the Google Groups
> "aroma.affymetrix" group.
> To unsubscribe from this group and stop receiving emails from it, send an
> email to aroma-affymetrix+unsubscr...@googlegroups.com.
> For more options, visit https://groups.google.com/groups/opt_out.
>
>
--
--
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/
---
You received this message because you are subscribed to the Google Groups
"aroma.affymetrix" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to aroma-affymetrix+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/groups/opt_out.
