Seqanswers Leaderboard Ad

Collapse

Announcement

Collapse
No announcement yet.
X
 
  • Filter
  • Time
  • Show
Clear All
new posts

  • #46
    If you're trying to do a comparison, the edgeR equivalent would be the "RLE" method in calcNormFactors(). estimateSizeFactors() doesn't take an option in this regard (see help(estimateSizeFactors) for what options it does have).

    Comment


    • #47
      Ah!
      Always a good answer from you. Thanks!

      Comment


      • #48
        Originally posted by Michael Love View Post
        The local fit helps you out by dipping down at the high counts, giving you lower estimates of dispersion there than by taking the mean overall. The large circles are genes whose dispersion estimates are not shrunk because they deviate too far from the fitted line. The local fit is useable as well I would say.

        I presume -- because you haven't posted the output of sessionInfo() -- that you could be using a version of DESeq2 which performs automatic independent filtering to optimize the number of rejections. Check the Details section of ?results and the argument 'independentFiltering'.
        Here is my sessionInfo()

        Code:
        > sessionInfo()
        R version 3.0.2 (2013-09-25)
        Platform: x86_64-apple-darwin10.8.0 (64-bit)
        
        locale:
        [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
        
        attached base packages:
        [1] grid      parallel  stats     graphics  grDevices utils     datasets  methods   base     
        
        other attached packages:
         [1] gplots_2.11.3           MASS_7.3-29             KernSmooth_2.23-10      caTools_1.14           
         [5] gdata_2.13.2            gtools_3.1.0            RColorBrewer_1.0-5      DESeq2_1.2.0           
         [9] RcppArmadillo_0.3.920.1 Rcpp_0.10.5             GenomicRanges_1.14.1    XVector_0.2.0          
        [13] IRanges_1.20.0          BiocGenerics_0.8.0     
        
        loaded via a namespace (and not attached):
         [1] annotate_1.40.0      AnnotationDbi_1.24.0 Biobase_2.22.0       bitops_1.0-6        
         [5] DBI_0.2-7            genefilter_1.44.0    lattice_0.20-24      locfit_1.5-9.1      
         [9] RSQLite_0.11.4       splines_3.0.2        stats4_3.0.2         survival_2.37-4     
        [13] tools_3.0.2          XML_3.95-0.2         xtable_1.7-1
        I'm still having some questions about the independent filtering. I am understanding that I have to perform it right after the:

        Code:
        dds<- estimateSizeFactors(dds)
        ddsLocal <- estimateDispersions(dds, fitType="local")
        ddsLocal <- nbinomWaldTest(ddsLocal)
        res <- results(ddsLocal, name= "condition_MR_vs_BR") 
        res <- res[order(res$padj),]
        With the code:

        Code:
        use <- res$baseMean >= 10 & !is.na(res$pvalue)
        resFilt <- res [use, ]
        resFilt$padj <- p.adjust(resFilt$pvalue, method="BH")
        And then continue with (for example):

        Code:
        #filter for upregulated and downregulated genes
        resSig <- resFilt[ resFilt$padj < .1, ]
        write.table(resSig[ order( resSig$log2FoldChange, -resSig$baseMean ), ]
                    ,"~\\DESeq2\\BRvsMR_old_DownRegulated.txt")
        write.table(resSig[ order( -resSig$log2FoldChange, -resSig$baseMean ),
                            ],"~\\DESeq2\\BRvsMR_old_UpRegulated.txt")
        But i also found this code:

        Code:
        filterThreshold <- 2.0
        keep <- rowMeans ( counts(ddsLocal, normalized=TRUE)) > filterThreshold
        And now i'm really confused. What is this for? (the last code) When i have to use it? Am I doing something wrong?

        Thank you again for the help.

        Comment


        • #49
          You are using DESeq2 v1.2 and reading the documentation off the Bioc website for version v1.0.

          This is why I recommended you to look up the man page for ?results.

          And it's always safer to use:

          browseVignettes("DESeq2")

          to ensure you are reading the documentation for the version you are using.

          Comment


          • #50
            Originally posted by Michael Love View Post
            You are using DESeq2 v1.2 and reading the documentation off the Bioc website for version v1.0.

            This is why I recommended you to look up the man page for ?results.

            And it's always safer to use:

            browseVignettes("DESeq2")

            to ensure you are reading the documentation for the version you are using.
            Ups! Thanks!

            By using:

            Code:
            browseVignettes("DESeq2")
            i am obtaining

            Code:
            /library/DESeq2/doc/DESeq2.pdf not found
            Any idea about why this is happening?

            Comment


            • #51
              Hi

              I just started using DEseq2, so please correct if i am wrong.
              My sample description has five time points (with 2 replicates).
              (t0,t0,t1,t1,t2,t2,t3,t3,t4,t4,t5,t5)

              if i set the timepoint design, i always get the differential expression through initial time point (control).

              In addition to this, i also want to have differential expression with respect to the immediate time points (like t1 to t0, t2 to t1, t3 to t2, t4 to t3 and t5 to t4). Can you please tell me how to set the factor levels to get this type of comparison.

              thank you
              Last edited by i4u412; 10-23-2013, 04:57 AM.

              Comment


              • #52
                Hi,

                I'm using DESeq2_1.2.0 and I am having this error:

                Code:
                > resSig <- res[ res$padj < .1, ]
                Error in normalizeSingleBracketSubscript(i, x, byrow = TRUE, exact = FALSE) : 
                  subscript contains NAs
                It's really weird because this error started today and I have been using the same datasets and the same Rscript for the last 3 weeks without this error coming up.

                Just in case, here is my code:

                Code:
                #Count matrix input
                Cele_SPvsLR_old = read.csv (file.choose(), header=TRUE, row.names=1)
                CeleDesign <- data.frame(
                  row.names = colnames(Cele_SPvsLR_old),
                  condition = factor(c("SP", "SP", "LR", "LR")))
                dds <- DESeqDataSetFromMatrix(countData = Cele_SPvsLR_old,
                                              colData = CeleDesign,
                                              design = ~ condition)
                dds
                
                #Est size factor = normalize for library size
                dds<- estimateSizeFactors(dds)
                ddsLocal <- estimateDispersions(dds, fitType="local")
                ddsLocal <- nbinomWaldTest(ddsLocal)
                plotDispEsts(ddsLocal)
                
                #Differential expression analysis
                resultsNames(ddsLocal)
                res <- results(ddsLocal, name= "condition_SP_vs_LR") 
                res <- res[order(res$padj),]
                head(res)
                plotMA(ddsLocal, ylim=c(-2,2), main="DESeq2")
                sum(res$padj < .1, na.rm=TRUE)
                
                #filter for significant genes 
                resSig <- res[ res$padj < 0.1, ]
                And my sessionInfo()

                Code:
                > sessionInfo()
                R version 3.0.2 (2013-09-25)
                Platform: x86_64-apple-darwin10.8.0 (64-bit)
                
                locale:
                [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
                
                attached base packages:
                [1] parallel  stats     graphics  grDevices utils     datasets  methods   base     
                
                other attached packages:
                [1] DESeq2_1.2.0            RcppArmadillo_0.3.920.1 Rcpp_0.10.5             GenomicRanges_1.14.1   
                [5] XVector_0.2.0           IRanges_1.20.0          BiocGenerics_0.8.0     
                
                loaded via a namespace (and not attached):
                 [1] annotate_1.40.0      AnnotationDbi_1.24.0 Biobase_2.22.0       DBI_0.2-7           
                 [5] genefilter_1.44.0    grid_3.0.2           lattice_0.20-24      locfit_1.5-9.1      
                 [9] RColorBrewer_1.0-5   RSQLite_0.11.4       splines_3.0.2        stats4_3.0.2        
                [13] survival_2.37-4      tools_3.0.2          XML_3.95-0.2         xtable_1.7-1
                Can anyone give me an explanation about why am I having this error now but not in the past using exactly the same Rscript and datasets? Does anyone know how to fix it?

                Thanks!
                Last edited by alisrpp; 10-24-2013, 11:50 AM.

                Comment


                • #53
                  Please use
                  Code:
                  resSig <- res[ which( res$padj < .1 ), ]
                  to get rid of the NAs.

                  Comment


                  • #54
                    Can I ask why you set it to .1 ?

                    Comment


                    • #55
                      contrast argument

                      Hi,

                      I have started using DESeq2 on a multi-factor design. Everything goes fine until I try to use the "contrast" argument on the "results" function in order to get somparisons other than the ones concerning the control condition.

                      This is what I get:
                      > res.R24.OE1=results(dds.gtype,contrast=c("gtype","R24","OE1"))
                      Error in results(dds.gtype, contrast = c("gtype", "R24", "OE1")) :
                      unused argument (contrast = c("gtype", "R24", "OE1"))

                      "gtype" being my factor and "R24" and "OE1" my two conditions.

                      I would be very grateful if you could help me out on this one.

                      Thanks!

                      Comment


                      • #56
                        Hi,
                        I am using DESeq2 version 1.2.5 and when I use this design:

                        > design(dds) <- formula(~ AGE + SEX + CONDITION)

                        Defining the variables:
                        > colData(dds)$CONDITION <- relevel(colData(dds)$CONDITION, "R")
                        > colData(dds)$SEX <- relevel(colData(dds)$SEX, "M")
                        Whereas AGE is a continuous variable

                        I get this error:
                        > dds <- DESeq(dds)
                        using pre-existing size factors
                        estimating dispersions
                        you had estimated dispersions, replacing these
                        gene-wise dispersion estimates
                        error: inv(): matrix appears to be singular

                        On the other hand, there are no problems with these combinations:
                        > design(dds) <- formula(~ AGE)
                        > design(dds) <- formula(~ CONDITION)
                        > design(dds) <- formula(~ SEX + CONDITION)

                        While I get the same error with:
                        > design(dds) <- formula(~ AGE + CONDITION)

                        Any clue what's going on here?
                        Thank you so much in advance for your help

                        Comment


                        • #57
                          Could "AGE" have been converted to a factor at some point?

                          Comment


                          • #58
                            Originally posted by dpryan View Post
                            Could "AGE" have been converted to a factor at some point?
                            No it has not. I checked:

                            > is(colData(dds)$AGE)
                            [1] "numeric" "vector" "atomic" "vectorORfactor"

                            Thanks for the suggestion

                            Comment


                            • #59
                              Can you please post colData(dds) (or some scrubbed version of it)? This helps me figure out what the design matrix looks like, and therefore why some designs are causing errors.

                              Comment


                              • #60
                                Originally posted by Michael Love View Post
                                Can you please post colData(dds) (or some scrubbed version of it)? This helps me figure out what the design matrix looks like, and therefore why some designs are causing errors.
                                Here it is:
                                > colData(dds)
                                DataFrame with 122 rows and 4 columns
                                CONDITION KEEP AGE SEX
                                <factor> <integer> <numeric> <factor>
                                ID1 A 1 47.06913 F
                                ID2 A 1 45.33333 M
                                ID3 A 1 59.63039 F
                                ID4 A 1 49.00753 M
                                ID5 A 1 34.94319 M
                                ... ... ... ... ...
                                ID118 B 1 30.00684 M
                                ID119 B 1 41.90828 F
                                ID120 B 1 39.04449 F
                                ID121 B 1 39.64682 M
                                ID122 B 1 16.70363 M

                                Comment

                                Latest Articles

                                Collapse

                                • seqadmin
                                  Choosing Between NGS and qPCR
                                  by seqadmin



                                  Next-generation sequencing (NGS) and quantitative polymerase chain reaction (qPCR) are essential techniques for investigating the genome, transcriptome, and epigenome. In many cases, choosing the appropriate technique is straightforward, but in others, it can be more challenging to determine the most effective option. A simple distinction is that smaller, more focused projects are typically better suited for qPCR, while larger, more complex datasets benefit from NGS. However,...
                                  10-18-2024, 07:11 AM
                                • seqadmin
                                  Non-Coding RNA Research and Technologies
                                  by seqadmin




                                  Non-coding RNAs (ncRNAs) do not code for proteins but play important roles in numerous cellular processes including gene silencing, developmental pathways, and more. There are numerous types including microRNA (miRNA), long ncRNA (lncRNA), circular RNA (circRNA), and more. In this article, we discuss innovative ncRNA research and explore recent technological advancements that improve the study of ncRNAs.

                                  Nobel Prize for MicroRNA Discovery
                                  This week,...
                                  10-07-2024, 08:07 AM

                                ad_right_rmr

                                Collapse

                                News

                                Collapse

                                Topics Statistics Last Post
                                Started by seqadmin, Yesterday, 05:31 AM
                                0 responses
                                10 views
                                0 likes
                                Last Post seqadmin  
                                Started by seqadmin, 10-24-2024, 06:58 AM
                                0 responses
                                20 views
                                0 likes
                                Last Post seqadmin  
                                Started by seqadmin, 10-23-2024, 08:43 AM
                                0 responses
                                48 views
                                0 likes
                                Last Post seqadmin  
                                Started by seqadmin, 10-17-2024, 07:29 AM
                                0 responses
                                58 views
                                0 likes
                                Last Post seqadmin  
                                Working...
                                X