Hot answers tagged

7

The prior count ends up getting scaled by the ratio of a library size to the average library size and then multiplied by 2 before getting added to each library size (I'm sure there's a good reason for that, but I don't know what it is). Using the example df data frame in your post, let's walk step-by-step through what cpm() is doing: df = data.frame(a = c(1,...


6

The more genes you have the more robust the scaling factor is (among the reason why one doesn't normalize to ERCC spike-ins without a compelling reason), so I suppose in theory it's better to filter after determining the scale factor. Having said that, I'd be surprised if the results changed much either way. Unless you end up filtering out a LOT of genes the ...


5

You are looking at the supplementary data of a paper. That seems to have given you a list of features, and some information about those features. Specifically, you seem to have a list of two types of element: "peg". Based on the information there, I assume "peg" stands for "protein encoding gene". Note that all peg lines have a ...


5

It is indeed true that for the DE analysis one should include batch into the formula to avoid changing the original counts. Still, for everything else such as plotting heatmaps use of removeBatchEffects is perfectly fine and (at least to me) a standard and well-accepted procedure. It essentially does not matter what you use to correct for the batch effect ...


4

In the DESeq2 manual there is a section titled "How can I include a continuous covariate in the design formula?" that deals with your question. Basically the process is no different from using a discrete covariate. DESeq <- DESeqDataSetFromMatrix(countData = counts, colData = metadata, design = ~ RIN + condition) In this example condition is the ...


4

FPKM are inherently experiment specific and can not be used to compare across samples. Let's consider the following two sequencing runs. Let $E1$ and $E2$ be the true, underlying expression in two samples of genes 1-6. Let $S1$ and $S2$ be the observed expression in our sequencing. $$ \begin{matrix} Gene & E1 & S1 & E2 & S2 \\ G1 & 100 &...


4

Not really an answer but an extended comment... and most likely something you don't like to hear I guess by technical replicates, it means taking the same biological sample and making 2 methylation libraries. If this is used as replicates in deseq2 or edger, the variance you are estimating is the technical variation that comes with preparing the library, ...


4

Ok, it turned out that some rows and/or columns contained all values as zeros and that was the reason for the error. I wrote a function to fix that: def fix_rawcounts(rawcounts, cell_ids, gene_symbols): '''The function checks whether there are any rows or columns with all zeros and removes them''' print('Checking whether rows or columns ...


3

You're really running out of degrees of freedom and the only actual replicate is a single sample from another patient, so please take the results with a huge grain of salt. You can use a design of ~cellType * coculture where cellType has levels NOF and CAF and coculture has levels yes and no (ideally you'd block by patient, but that's not an option for you)....


3

In general, you would expect roughly the same distribution on both sides of your volcano plot. The first two plots are concerning since there's so much on the right and not on the left. With such a small number of samples, it is hard to conclude if this is an artifact of not having enough samples, a technical artifact from the way the data was processed, or ...


3

Random partitioning (as you have done) seems like a reasonable thing to do to work out what's going on. From the results you have suggested, it looks like DESeq2 is performing better than EdgeR, so use DESeq2. Can you show the volcano plots and the partition sizes? Note: I'm a little bit biased in this, because DESeq2 is what I prefer to use (particularly ...


3

The help page for removeBatchEffect() explicitly states that it should not be used in conjunction with differential expression testing, so don't do that. You already have the batch effect included in your model, so you're already correcting for it when testing for differential expression. If your MDS plots suggest that samples are variably affected by the ...


3

If you want the values by group then subset y to contain the samples of interest and then feed that to aveLogCPM().


2

Samples will only cluster by experimental group if the experimental effect is large enough that it's the primary source of variance between your samples. If that's not the case then you'll get results like you're seeing. Whether that's a problem or not ends up depending on how large of an effect you expect to see. If your groups are different stages in ...


2

I can see that there would be a three step process to doing this: Merge counts from all samples in the group and then resample pseudo-replicates from this. If x is a matrix with samples in the group being columns and genes being rows s <- sample(row.names(x), n = mean(colSums(x)), probs=rowSums(x)/sum(rowSums(x)) stab <- table(s) s <- as.vector(...


2

If I'm interested to get normalized gene counts, can I go ahead and multiple each individual's norm.factor into its gene counts? For example, the expected count for IND1 for Gene1 is 100 and it's norm.factor is 0.80, can I say that the normalized gene count is 100*0.80=80? No; you divide, but if you poke around, you can probably find a way to get ...


2

The simplest thing to do is to trim the 150 bp fasts so that they are 100 long. I don't think there is an easy way to correct for the fact that the 150 bp long reads will have a higher unambiguous rate of alignment and gene assignment than the 100 bp long reads. If you have a mix of all the experimental conditions across all the lengths, you can include ...


1

It makes no difference if you process the BAM files one at a time with featureCounts or all together, except that it changes how you have to read the files into R. You can supply edgeR with lists of contrasts to have it compute fold-changes and p-values for. Please have a look at the edgeR user guide for examples.


1

If you got DEGs with statistically significant p-values after multiple testing following a valid pipeline, I'd believe that, even if looking at the first few PCs doesn't look promising.


1

You do the normalization before running your edgeR. The purpose of RUVg is to remove "Remove Unwanted Variation Using Control Genes". In your code, you ran edgeR and then normalize the data using RUVg, which is only going to return you the normalized counts. Using the example dataset in vignette: library(RUVSeq) library(zebrafishRNASeq) data(...


1

You may run edgeR for methylation analysis without replicates (https://f1000research.com/articles/6-2055). But I also recommend you to have a look at the R DSS package (http://bioconductor.org/packages/release/bioc/html/DSS.html). It has a smoothing step which allows to modelize the biological variability by taking information from neighboring sites.


1

gi looks like a Genbank Identifier. You can search for the 10 digit number following gi| at the Genbank website.


1

I would use the p-value/FDR which each method returns to rank the gene list for that method in order from 'most likely to be DE' to least likely. You would then have three ranked lists of genes - you would put the actual rankings in for each gene - 1 for the top ranked, 2 for the next, etc. You could then carry out, say Spearman's rank correlation on the ...


1

You'll have a bit more power if you use all of the samples at once. Nothing is useful for this, don't waste your time with it. You're not looking for outliers, you're looking to see how samples cluster. You can just use correlation for that.


1

Yes this is nothing but a z-score. The scale function calculates column z-scores, but since you use t (transpose), it is actually performed on rows. There are many other ways to calculate z-scores, but if done correctly they will all give the same results.


Only top voted, non community-wiki answers of a minimum length are eligible