Hot answers tagged

8

If you compare A vs B the genes's fold change will have the opposite sign to B vs A. So will be the gene set up or down-regulated depending on the comparison to take. The gene set test analyze if a given group of elements is sorted in a certain way in a list (I am talking about the GSE like the one you performed or the one in the Broad Institute). Usually ...


5

Seruat will give you a list of genes which it thinks are upregulated in a particular cluster. Look at the functions that talk about marker genes - these functions basically do a DE analysis of the genes in one cluster compared to the others. Then take that list and feed it to any standard GO analysis tool. Have a look at the topGO topKEGG and geneSetTest ...


5

While I'm not familiar with GSEA software in particular, I believe your problem is that it only tests for upregulated gene sets. Notice: 866/4408 gene sets are upregulated 3542/4408 gene sets are upregulated 3542+866 = 4408. I.e., 866 sets have higher mean expression in positive condition, the rest have higher mean expression in negative ...


4

I assume that you are talking about the implementation of these methods in the limma package. Otherwise this answer does not apply. I think that your questions can be answered with some simulations where we can test with some "genes" with a known relationship: library("limma") set.seed(123) # Create some genes and samples nGenes <- 40 nSamples <- ...


3

I have answered my own question however, it is specific to windows users.. If you get the same error from the cmd as I did, you must manually put in environment variables related to Java. Do as the user Mohamed did in this SO post: https://stackoverflow.com/questions/9303889/error-occurred-during-initialization-of-vm-could-not-reserve-enough-space-for ...


3

If using DESeq2 with GSEA, I'd recommend ranking by shrunk log2FC values. It'd also be worth considering ranking positive and negative associations separately, because the standard GSEA algorithm doesn't cross at the zero point when associations change from positive to negative. You shouldn't be using p-values to rank anything. A p-value gives information ...


3

Ok. Answering my own question. I did this way. And I got the output I want. install.packages("msigdbr") library(msigdbr) m_df = msigdbr(species = "Homo sapiens", category = "C5") head(m_df) Output: gs_name gs_id gs_cat gs_subcat human_gene_symb~ species_name entrez_gene <chr> <chr> <chr> <chr> <chr> <...


3

I finally found another answer to my question. Please read this great article in May 12 2017 BioMed Central (BMC) Bioinformatics article titled Ranking metrics in gene set enrichment analysis: do they matter?. Also, please read this blog , Diving into Genetics and Genomics: Gene Set Enrichment Analysis (GSEA) explained. After reading these two articles, my ...


3

You seem to have forgotten the space after the #. It should be as in here (Just for future readers know where to find the specs). Also it should be a single line (Maybe this error was introduced while pasting here or it has been a duplicate because there are 39 samples in the 3 line and they should be 37). Either adding 2 more samples to the initial count ...


2

Have you checked firewall settings on your laptop? See Firewall / FTP connection issues at: GSEA wiki - known issues


2

I tinkered with a program GSEA-SNP quite a few years back, which claims that it does a similar ranking procedure with SNPs. It carries out its procedure by first linking SNPs with genes, then running an algorithm similar to GSEA. There's a bit more detail in the paper. Unfortunately this is a space I've got a bit further away from in recent years, so I don'...


2

You rank the fit coefficient rather than the original score matrix. So, given a score matrix, D: D = matrix(c(22,20,9,8,46,22,18,10,3,18,3,29,2,1,5,45,43,47,17,5,14,44,21,36), byrow=T, ncol=6) cl = c(0,0,0,1,1,1) s = apply(m, 1, function(x) coef(lm(x~cl))[2]) # [1] o = rank(s) o = max(o) - o + 1 # [2] o is then the rank of each row. [1] This fits each row ...


2

There is no purpose-built R package to perform gene set enrichment analysis on single-cell data but there does not need to be. You should be able to tools developed for bulk-RNA-Seq or microarray data, although you may not get as significant results from a sparse scRNA-Seq matrix as single-cell technologies have poor sensitivity and miss genes. What you need ...


2

I personally rank by signed nominal p-value. The fgsea author recommended on biostars.org to use either -log10(nominal p-value) or the F-statistics column (or whatever statistic the tool you use outputs) but not FDR/adjusted p-values as the latter produces a lot of ties for genes with low significances.


2

Each of these methods do something different, so it is reasonable to expect different results. The bottom line is that there isn't a single question you can ask when you do an "enrichment test". For example: you can test if a set of genes are more expressed than this other list of genes or they are above the mean levels of any other gene. Even more,...


2

You could apply the GSVA transformation to later compare them. You would need to give a step back and start from the expression values not log2FC, but I think it would be better than using the log2FC directly. If you need to use them you could compare the gene set enrichment value of these two sets on the same comparison. I would use fgsea with the log2FC ...


2

I'm not sure I understand "I would like to run GSEA or a similar analysis to find WGCNCA clusters that differ based on the interaction between the two main variables". I would run an association analysis (regression) for module eigengenes and the appropriate interaction term. The (most) significant modules are your candidates. This step is simply ...


2

Are you calling java with sufficient memory? Try increasing -Xmx2G to -Xmx8G, provided your PC has 8 GB RAM.


2

Actually I believe I get it now. i is not a constant. We vary it from i=1 to N and consider the i for which $P_{hit} - P_{miss}$ is largest in magnitude


1

I'd probably just use fgsea. Your code won't have to change much, it also just uses a ranked list.


1

The noise is comming from the subset.mask, which is created above in a loop with the number of permutations. for (r in 1:nperm) { #L90 .... subset.mask[, r] <- as.numeric(c(subset.class1, subset.class2)) # L107 So by multiplying the random selection of subsets by the expression we get the "noise". Later on line 251 we get the ratio of noise/...


1

In Deseq2, it depends on contrast eg. if your contrast <- c("Condition", firstC, SecondC) is then -ve is downregulated in FirstC and same for +ve. It is always better to take look at count data for the confirmation (responders versus non-responders samples). [for more detail you can refer here][1] [1]https://github.com/amarinderthind/RNA-seq-...


1

You should be using the most positive fold changes, not the most negative. The positive ones are the ones that are most highly expressed in your cluster of interest. You probably want to use less than 100 genes. Depending on the experiment, you may not have that many good markers for a particular cluster. The other question you mentioned refers to DESeq2 ...


1

Yes, if you have a ranked list, you can use GSEA to calculate the enrichment of the strains along the ranked list. Note that hypergeometric test is for count data, so you would test if there are more strains of one type that fit on certain conditions. I would recommend this article comparing different enrichment methods to understand well what they do.


1

You can use GSA.read.gmt function from GSA package. The following code can be used to convert the file to a dataframe. Just ignore the warnings. Original_response library(GSA) data <- GSA.read.gmt("c5.all.v6.2.symbols.gmt") gene_names <- unlist(data$genesets, use.names=FALSE) your_dataframe <- cbind(data$geneset.names,gene_names) colnames(...


1

You can import the csv file as a table. This does not require further libraries: #you're reading a csv file, using tab as field separator and considering the first line as the headers of the data table data <- read.csv(filename, sep="\t", header=TRUE) then what you need are just two columns of your data: colnames(data) #will give you the names of the ...


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