pachterlab / sleuth

Differential analysis of RNA-Seq
http://pachterlab.github.io/sleuth
GNU General Public License v3.0
305 stars 95 forks source link

Trouble with the Sleuth_wt function #154

Closed ScottSchumacker closed 6 years ago

ScottSchumacker commented 6 years ago

Hello! I am currently trying to test for differentially expressed genes with Sleuth. I am able to create a likely hood ratio test but am struggling with creating a wald test.

Here is how I set up the models for my data: so <- sleuth_fit(so, ~tissue, 'full')
so <- sleuth_fit(so, ~1, 'reduced')

When I use the models(so) function, I get the following output:

 formula:  ~tissue 
 coefficients:
     (Intercept)
     tissueretina
 [  reduced  ]
 formula:  ~1 
 coefficients:
     (Intercept)

I tried using the sleuth_wt function to create a wald test and here is the error I get: sleuth_wt(so, which_beta = 'tissueretina', which_model = 'full')

     sleuth object

 bears: 18 
 design: NULL

Any help is appreciated! Thank you!

pimentel commented 6 years ago

Hi @BIOProgrammer21

Could you please post the full error message from sleuth_wt? Everything you are doing looks correct. Also, which version of sleuth are you running?

Thanks,

Harold

pimentel commented 6 years ago

Actually, to be clear, could you post your entire sleuth_wt call?

It should look like this:

so <- sleuth_wt(so, which_beta = 'tissueretina', which_model = 'full')
the_results <- sleuth_results(so, 'tissueretina', 'wt')

Best,

Harold

ScottSchumacker commented 6 years ago

Hi @pimentel,

Thank you for the response. Here is the output from the sleuth_wt call as well as the output for calling models(so) and tests(so).

 > sleuth_wt(so, which_beta = 'tissueretina', which_model = 'full')
     sleuth object

 bears: 18 
 design: NULL 

 > the_results <- sleuth_results(so, 'tissueretina', 'wt')
 Error in get_test(obj, test, "wt", which_model) : 
   'tissueretina' is not a valid label for a test. Please see valid models and tests using the functions       'models' and 'tests'. Remember to also correctly specify the test type.

 > models(so)
 [  full  ]
 formula:  ~tissue 
 coefficients:
     (Intercept)
     tissueretina
 [  reduced  ]
 formula:  ~1 
 coefficients:
     (Intercept)

 > tests(so)
 ~likelihood ratio tests:
     reduced:full

 ~wald tests:
     no tests found.
 > 

As for the version of Sleuth I am using, I am not sure. I followed the instructions here https://pachterlab.github.io/sleuth/download when downloading it

warrenmcg commented 6 years ago

I have identified the problem, and now I understand the output from sleuth_wt. You'll notice that you when run tests(so), it doesn't show any Wald tests. You need to assign the results from sleuth_wt to so, like the following:

so <- sleuth_wt(so, which_beta = 'tissueretina', which_model = 'full')

Then you can run the_results <- sleuth_results(so, 'tissueretina', 'wt').

The important part is so <-. Without it, R doesn't know that you updated so to include the Wald test. The output below the function call is R printing out a summary of the object that should have been assigned:

    sleuth object

 bears: 18 
 design: NULL 

There are some functions that would modify the so object in place (without assignment), but those are unusual in R. Usually, if you need to update the object, you assign the function results to the object. That's how all of sleuth's functions work.

A simple example is to make a data frame versus assigning it to an object:

data.frame(numbers = 1:10, letters = rep("A", 10)) # this will print the result to the command line
df <- data.frame(numbers = 1:10, letters = rep("A", 10)) # this will store the result to the object 'df'

Hope that helps! Because we figured out the issue, and it isn't a bug, I'm closing this issue.

ScottSchumacker commented 6 years ago

@warrenmcg Thank you! I am surprised I missed that.

I appreciate the help!

warrenmcg commented 6 years ago

@BIOProgrammer21 you and me both.