git-svn-id: file:///home/git/hedgehog.fhcrc.org/bioconductor/trunk/madman/Rpacks/crlmm@48955 bc3139a8-67e5-0310-9ffc-ced21a209358
... | ... |
@@ -769,8 +769,8 @@ fit.lm1 <- function(idxBatch, |
769 | 769 |
marker.index, |
770 | 770 |
object, |
771 | 771 |
Ns, |
772 |
- normal, |
|
773 |
- snpflags, |
|
772 |
+## normal, |
|
773 |
+## snpflags, |
|
774 | 774 |
batchSize, |
775 | 775 |
SNRMin, |
776 | 776 |
MIN.SAMPLES, |
... | ... |
@@ -977,8 +977,8 @@ fit.lm2 <- function(idxBatch, |
977 | 977 |
marker.index, |
978 | 978 |
object, |
979 | 979 |
Ns, |
980 |
- normal, |
|
981 |
- snpflags, |
|
980 |
+## normal, |
|
981 |
+## snpflags, |
|
982 | 982 |
batchSize, |
983 | 983 |
SNRMin, |
984 | 984 |
MIN.SAMPLES, |
... | ... |
@@ -1088,8 +1088,8 @@ fit.lm3 <- function(idxBatch, |
1088 | 1088 |
marker.index, |
1089 | 1089 |
object, |
1090 | 1090 |
Ns, |
1091 |
- normal, |
|
1092 |
- snpflags, |
|
1091 |
+## normal, |
|
1092 |
+## snpflags, |
|
1093 | 1093 |
batchSize, |
1094 | 1094 |
SNRMin, |
1095 | 1095 |
MIN.SAMPLES, |
... | ... |
@@ -1327,8 +1327,8 @@ fit.lm4 <- function(idxBatch, |
1327 | 1327 |
marker.index, |
1328 | 1328 |
object, |
1329 | 1329 |
Ns, |
1330 |
- normal, |
|
1331 |
- snpflags, |
|
1330 |
+## normal, |
|
1331 |
+## snpflags, |
|
1332 | 1332 |
batchSize, |
1333 | 1333 |
SNRMin, |
1334 | 1334 |
MIN.SAMPLES, |
... | ... |
@@ -1442,7 +1442,6 @@ fit.lm4 <- function(idxBatch, |
1442 | 1442 |
if(any(pseudoAR)){ |
1443 | 1443 |
nu1[pseudoAR] <- 2^(mus[pseudoAR, 1]) - 2*phi1[pseudoAR] |
1444 | 1444 |
} |
1445 |
- |
|
1446 | 1445 |
## normal.f <- NORM.np[, k] |
1447 | 1446 |
## A.F <- A.F*normal.f[, gend==2] |
1448 | 1447 |
A.F[A.F==0] <- NA |
... | ... |
@@ -2774,6 +2773,7 @@ ellipseCenters <- function(object, index, allele, batch, log.it=TRUE){ |
2774 | 2773 |
return(centers) |
2775 | 2774 |
} |
2776 | 2775 |
|
2776 |
+ |
|
2777 | 2777 |
computeCN <- function(object, |
2778 | 2778 |
filenames, |
2779 | 2779 |
which.batches, |
... | ... |
@@ -2798,89 +2798,44 @@ computeCN <- function(object, |
2798 | 2798 |
stopifnot("position" %in% fvarLabels(object)) |
2799 | 2799 |
stopifnot("isSnp" %in% fvarLabels(object)) |
2800 | 2800 |
batch <- batch(object) |
2801 |
+ is.snp <- isSnp(object) |
|
2802 |
+ is.autosome <- chromosome(object) < 23 |
|
2803 |
+ is.annotated <- !is.na(chromosome(object)) |
|
2804 |
+ is.X <- chromosome(object) == 23 |
|
2801 | 2805 |
if(type=="autosome.snps"){ |
2802 |
- marker.index <- which(chromosome(object) < 23 & isSnp(object) & !is.na(chromosome(object))) |
|
2806 |
+ marker.index <- which(is.snp & is.autosome & is.annotated) |
|
2803 | 2807 |
nr <- length(marker.index) |
2804 | 2808 |
Ns <- initializeBigMatrix("Ns", nr, 5) |
2805 | 2809 |
colnames(Ns) <- c("A", "B", "AA", "AB", "BB") |
2806 |
- if(!file.exists(file.path(ldPath(), "normal.snps.rda"))){ |
|
2807 |
- normal <- initializeBigMatrix("normal.snps", nr, ncol(object), vmode="integer") |
|
2808 |
- normal[,] <- 1L |
|
2809 |
- save(normal, file=file.path(ldPath(), "normal.snps.rda")) |
|
2810 |
- } else load(file.path(ldPath(), "normal.snps.rda")) |
|
2811 |
- if(!file.exists(file.path(ldPath(), "flags.snps.rda"))){ |
|
2812 |
- flags <- initializeBigMatrix("flags.snps", nr, length(unique(batch(object))), vmode="integer") |
|
2813 |
- flags[,] <- 0L |
|
2814 |
- save(flags, file=file.path(ldPath(), "flags.snps.rda")) |
|
2815 |
- } else{ |
|
2816 |
- load(file.path(ldPath(), "flags.snps.rda")) |
|
2817 |
- } |
|
2818 | 2810 |
if(verbose) message("Estimating allele-specific copy number at autosomal SNPs") |
2819 | 2811 |
FUN <- fit.lm1 |
2820 | 2812 |
} |
2821 | 2813 |
if(type=="autosome.nps"){ |
2822 |
- marker.index <- which(chromosome(object) < 23 & !isSnp(object) & !is.na(chromosome(object))) |
|
2823 |
- ##marker.index <- which(chromosome(object) < 23 & isSnp(object) & !is.na(chromosome(object))) |
|
2814 |
+ marker.index <- which(is.autosome & !is.snp & is.annotated) |
|
2824 | 2815 |
nr <- length(marker.index) |
2825 | 2816 |
Ns <- initializeBigMatrix("Ns", nr, 5) |
2826 | 2817 |
colnames(Ns) <- c("A", "B", "AA", "AB", "BB") |
2827 |
- if(!file.exists(file.path(ldPath(), "normal.snps.rda"))){ |
|
2828 |
- normal <- initializeBigMatrix("normal.nps", nr, ncol(object), vmode="integer") |
|
2829 |
- normal[,] <- 1L |
|
2830 |
- save(normal, file=file.path(ldPath(), "normal.nps.rda")) |
|
2831 |
- } else load(file.path(ldPath(), "normal.nps.rda")) |
|
2832 |
- if(!file.exists(file.path(ldPath(), "flags.nps.rda"))){ |
|
2833 |
- flags <- initializeBigMatrix("flags.nps", nr, length(unique(batch(object))), vmode="integer") |
|
2834 |
- flags[,] <- 0L |
|
2835 |
- save(flags, file=file.path(ldPath(), "flags.nps.rda")) |
|
2836 |
- } else{ |
|
2837 |
- load(file.path(ldPath(), "flags.nps.rda")) |
|
2838 |
- } |
|
2839 | 2818 |
if(verbose) message("Estimating allele-specific copy number at autosomal SNPs") |
2840 | 2819 |
FUN <- fit.lm2 |
2841 | 2820 |
} |
2842 | 2821 |
if(type=="X.snps"){ |
2843 |
- marker.index <- which(chromosome(object) == 23 & isSnp(object) & !is.na(chromosome(object))) |
|
2844 |
- ##marker.index <- which(chromosome(object) < 23 & isSnp(object) & !is.na(chromosome(object))) |
|
2822 |
+ marker.index <- which(is.X & is.snp) |
|
2845 | 2823 |
nr <- length(marker.index) |
2846 | 2824 |
Ns <- initializeBigMatrix("Ns", nr, 5) |
2847 | 2825 |
colnames(Ns) <- c("A", "B", "AA", "AB", "BB") |
2848 |
- if(!file.exists(file.path(ldPath(), "normal.X.snps.rda"))){ |
|
2849 |
- normal <- initializeBigMatrix("normal.X.snps", nr, ncol(object), vmode="integer") |
|
2850 |
- normal[,] <- 1L |
|
2851 |
- save(normal, file=file.path(ldPath(), "normal.X.snps.rda")) |
|
2852 |
- } else load(file.path(ldPath(), "normal.X.snps.rda")) |
|
2853 |
- if(!file.exists(file.path(ldPath(), "flags.X.snps.rda"))){ |
|
2854 |
- flags <- initializeBigMatrix("flags.X.snps", nr, length(unique(batch(object))), vmode="integer") |
|
2855 |
- flags[,] <- 0L |
|
2856 |
- save(flags, file=file.path(ldPath(), "flags.X.snps.rda")) |
|
2857 |
- } else{ |
|
2858 |
- load(file.path(ldPath(), "flags.X.snps.rda")) |
|
2859 |
- } |
|
2860 | 2826 |
if(verbose) message("Estimating allele-specific copy number at autosomal SNPs") |
2861 | 2827 |
FUN <- fit.lm3 |
2862 | 2828 |
} |
2863 | 2829 |
if(type=="X.nps"){ |
2864 |
- marker.index <- which(chromosome(object) == 23 & !isSnp(object) & !is.na(chromosome(object))) |
|
2830 |
+ marker.index <- which(is.X & !is.snp) |
|
2865 | 2831 |
nr <- length(marker.index) |
2866 | 2832 |
Ns <- initializeBigMatrix("Ns", nr, 5) |
2867 | 2833 |
colnames(Ns) <- c("A", "B", "AA", "AB", "BB") |
2868 |
- if(!file.exists(file.path(ldPath(), "normal.X.nps.rda"))){ |
|
2869 |
- normal <- initializeBigMatrix("normal.X.nps", nr, ncol(object), vmode="integer") |
|
2870 |
- normal[,] <- 1L |
|
2871 |
- save(normal, file=file.path(ldPath(), "normal.X.nps.rda")) |
|
2872 |
- } else load(file.path(ldPath(), "normal.X.nps.rda")) |
|
2873 |
- if(!file.exists(file.path(ldPath(), "flags.X.nps.rda"))){ |
|
2874 |
- flags <- initializeBigMatrix("flags.X.nps", nr, length(unique(batch(object))), vmode="integer") |
|
2875 |
- flags[,] <- 0L |
|
2876 |
- save(flags, file=file.path(ldPath(), "flags.X.nps.rda")) |
|
2877 |
- } else{ |
|
2878 |
- load(file.path(ldPath(), "flags.X.nps.rda")) |
|
2879 |
- } |
|
2880 | 2834 |
if(verbose) message("Estimating allele-specific copy number at autosomal SNPs") |
2881 | 2835 |
FUN <- fit.lm4 |
2882 | 2836 |
} |
2883 | 2837 |
index.strata <- splitIndicesByLength(marker.index, ocProbesets()) |
2838 |
+ ## create a separate object for each strata and combine() at the very end. |
|
2884 | 2839 |
obj <- construct(filenames=filenames, |
2885 | 2840 |
cdfName=annotation(object), |
2886 | 2841 |
copynumber=TRUE, |
... | ... |
@@ -2892,8 +2847,8 @@ computeCN <- function(object, |
2892 | 2847 |
marker.index=marker.index, |
2893 | 2848 |
object=obj, |
2894 | 2849 |
Ns=Ns, |
2895 |
- normal=normal, |
|
2896 |
- snpflags=flags, |
|
2850 |
+ ##normal=normal, |
|
2851 |
+ ##snpflags=flags, |
|
2897 | 2852 |
snpBatches=index.strata, |
2898 | 2853 |
batchSize=ocProbesets(), |
2899 | 2854 |
SNRMin=SNRMin, |