Dear Henrik, thans for the useful answers!
RAS = Relative Allele Scorel is normally used as a frequency estimate for pooled DNA experiments. I saw around that SNP6.0 seems not to be used for this scope before, so there are not really validated methods to estimate frequencies from it. 1. The RAS score is defined as RAS = I_B/(I_A+I_B) , so it seems to me to be equal to "fracB" as long as I am correct to assume that for the snps, using combineAlleles=FALSE, thetaA = median(PMA) and thetaB=median(PMB). This is my understanding from the CRMA2 paper. Using the notation of the paper: thataA_i_j = median_k(Y_i_j_k_A) 2. By allelc crosstalk normalization I obtain quite a number of negative probe intensities (around 3%). Do you have nay suggestion about how to handle those? I guess this is a general problem since from my understanding the matrix S is calculated without constraints. 3. Is there a way to upload figures on this group? I have done some diagnositc plots that might be interesting for other users. Best Regards Marco On 26 Juni, 07:35, Henrik Bengtsson <[email protected]> wrote: > Hi. > > 2009/6/24 marco <[email protected]>: > > > > > Dear Henrik, > > > I wonder if aroma.affymetrix is suitable for pooledGWAstudies and > > if you have any experience with this subject. In principle this > > consists in extracting the RAS scores, so should be doable. I have the > > following specific questions and snippet: > > > 1. by using : > > > plm <- AvgCnPlm(csN, mergeStrands=TRUE, combineAlleles=FALSE) > > data <- extractTotalAndFreqB(ce, units=units); > > > - the freqB column should be the equivalent of RAS? > > Yes, there should be a one-to-one relationship between fracB and RAS. > I don't have the definition of RAS in front of me, but fracB is > defined as fracB = thetaB / (thetaA+thetaB) [non logged signals!] so > you should be able to find the transform. > > (Nowadays, we prefer to call it allele-B fraction instead of allele-B > frequency; the package will later reflect this too). > > > - how is it possible to get the SNPs names in the data.frame from > > extractTotalAndFreqB > > Note, that function does *not* return a 'data.frame' but a 'matrix' or > an 'array' (see below). > > unitNames <- getUnitNames(cdf, units=units); > > > - can the freqB for all the arrays be extracted at once, or only array > > by array? > > fracB <- extractTotalAndFreqB(ces, units=units); > > returns a Jx2xI array where J=length(units) and I=nbrOfArrays(ces), whereas, > > ce <- getFile(ces, 1); > fracB <- extractTotalAndFreqB(ce, units=units); > > returns a Jx2 matrix. > > > > > 2. Is there a way to eliminate CN probes from the analysis from the > > very beginning (just SNPs are need for pooling)? > > Don't know what you mean from the very beginning, but you can always > find the types of a CDF unit using: > > unitTypes <- getUnitTypes(cdf); > > SNPs have value 2 (this is according to the Affymetrix CDF file format), e.g. > > snpUnits <- which(unitTypes == 2); > > will identify the SNP units for you. > > > > > 3. Will be the standard normalization steps of CRMA usefuls/ > > deleterious/neutral? > > Can you restate this question; not sure I understand it. > > > > > I eliminated the crosstalk step so far, since I am not sure how this > > behaves on a DNA pool. It seems to me that the crosstalk requires 3 > > clusters for AA, AB BB to work properly. Is this correct? > > Please see CRMA and CRMA v2 papers [1,2] which explains this. The > crosstalk estimate should be fairly robust against what is going on > inside the genotype "cone". It tries to fit the crosstalk model by > giving much more weight to data points there the (A,AA,AAA,...) and > (B,BB,BBB,...) radials. > > > > > 4. Same question for the pcr and probe sequence normalization steps. > > Seems to me that they should not disturb > > They should be applicable as is. > > > > > > > 5. Do you see any missing/unuseful/wrong step in the code below? > > > ############# > > cdf <- AffymetrixCdfFile$byChipType("GenomeWideSNP_6", > > tags="Full") > > gi <- getGenomeInformation(cdf) > > si <- getSnpInformation(cdf) > > acs <- AromaCellSequenceFile$byChipType(getChipType(cdf, > > fullname=FALSE)) > > cdf <- AffymetrixCdfFile$byChipType("GenomeWideSNP_6", > > tags="Full") > > csR <- AffymetrixCelSet$byName("POOLS_SNP", cdf=cdf) > > > #acc <- AllelicCrosstalkCalibration(csR, model="CRMAv2") > > #print(acc) > > #csC <- process(acc, verbose=verbose) > > #print(csC) > > > bpn <- BasePositionNormalization(csR, target="zero") > > csN <- process(bpn, verbose=verbose) > > plm <- AvgCnPlm(csN, mergeStrands=TRUE, combineAlleles=FALSE) > > > if (length(findUnitsTodo(plm)) > 0) { > > units <- fitCnProbes(plm, verbose=verbose) > > str(units) > > units <- fit(plm, verbose=verbose) > > str(units) > > } > > ces <- getChipEffectSet(plm) > > fln <- FragmentLengthNormalization(ces, target="zero") > > cesN <- process(fln, verbose=verbose) > > > ce <- getFile(cesN, 1) > > data <- extractTotalAndFreqB(ce, units=units); > > Look alright to me. > > REFERENCES: > [1] H. Bengtsson; R. Irizarry; B. Carvalho; T.P. Speed, Estimation and > assessment of raw copy numbers at the single locus level, > Bioinformatics, 2008. [pmid: 18204055] > [2] H. Bengtsson; P. Wirapati & T.P, Speed, A single-array > preprocessing method for estimating full-resolution raw copy numbers > from all Affymetrix genotyping arrays including GenomeWideSNP 5 & 6, > Biinformatics, 2009. > > /Henrik > > PS. Please try to split up unrelated questions in separate messages. > > --~--~---------~--~----~------------~-------~--~----~ 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. To post to this group, send email to [email protected] To unsubscribe from this group, send email to [email protected] For more options, visit this group at http://groups.google.com/group/aroma-affymetrix?hl=en -~----------~----~----~----~------~----~------~--~---
