Methylkit usage for 'unite'

[Kalyan K Pasumarthy]

Hi,

I am using methylkit for RRBS data. I have 19 samples as test/treatment_0 and 15 samples as control/treatment_1. I want to perform 'unite' function on these samples. However, the way i want to use this function is different from the way it works.

I want to apply the unite function on equal percentage of samples. For example, I want CpGs which are represented in 80% of the samples (=15 in test samples and 10 in control, numbers were rounded off). However, to my understanding methylkit only allows 'unite' command to function in terms of numbers equally applicable on both treatment conditions ie., I have only option to choose  a definite number of samples which contain a particular CpG (ie., I have to choose min.pergroup=14L and this implies 73% samples in test condition and 93% in control)

Is there a way where I can choose to run 'unite' function to be applicable on equal percentage of samples(ie., min.pergroup = 80 which means 15 samples in test and 10 samples in control) ? In the other case can I choose to apply 'unite' function with separate 'min.pergroup' selection on each treatment (ie., min.pergroup.treatment0=15 and min.pergroup.treatment1=10).

Regards,
Kalyan

[Altuna Akalin]

there is no way to do that right now, and at the moment we are not planning to implement this functionality. But it is not hard to code, if you have experience in R. 

Here are couple of alternatives how to achieve your goal:

1)You can check out the code in github, add this functionality and make a pull request or just use it for your own purpose. 

2)Alternatively, you can just write a new unite() function let's call it unite2() that implements this functionality as an R function that is not part of the package. Example of such a function is here (just an example, it is not doing what you need to do):

https://gist.github.com/al2na/a9ea925fc11af26b9c72

3) use the unite function in the package with min.per.group=10L, then filter the methylBase object on the test group of samples so that the test group has at least 15 samples covered for a given CpG. 

Best,

Altuna

--
You received this message because you are subscribed to the Google Groups "methylkit_discussion" group.
To unsubscribe from this group and stop receiving emails from it, send an email to methylkit_discussion+unsub...@googlegroups.com.
To post to this group, send email to methylkit_discussion@googlegroups.com.
Visit this group at http://groups.google.com/group/methylkit_discussion.
For more options, visit https://groups.google.com/groups/opt_out.

[Kalyan K Pasumarthy]

Hi Altuna,

Thanks for letting me know the multiple options. My familiarity with R (and git) is limited to basics and I can understand the packages and use the functions.

I was trying to figure out the location of code for 'unite' function on github. While I was able to see the code for other functions at https://github.com/al2na/methylKit/tree/master/R, I could not locate the code for 'unite'.

I went through the newUnite() function and I could understand the source.

I would appreciate if you could point me to exact source for 'unite' function from methylkit(0.9.2), so that I can modify it to suit me!

Regards,
Kalyan

On 7 February 2014 14:24, Altuna Akalin <aak...@gmail.com> wrote:

there is no way to do that right now, and at the moment we are not planning to implement this functionality. But it is not hard to code, if you have experience in R. 

Here are couple of alternatives how to achieve your goal:

1)You can check out the code in github, add this functionality and make a pull request or just use it for your own purpose. 

2)Alternatively, you can just write a new unite() function let's call it unite2() that implements this functionality as an R function that is not part of the package. Example of such a function is here (just an example, it is not doing what you need to do):

https://gist.github.com/al2na/a9ea925fc11af26b9c72

3) use the unite function in the package with min.per.group=10L, then filter the methylBase object on test group to make so that the test group has at least 15 samples covered for a CpG. 

Best,

Altuna

On Fri, Feb 7, 2014 at 12:39 PM, Kalyan K Pasumarthy <kaly...@gmail.com> wrote:

--

You received this message because you are subscribed to the Google Groups "methylkit_discussion" group.

To unsubscribe from this group and stop receiving emails from it, send an email to methylkit_discus...@googlegroups.com.
To post to this group, send email to methylkit_...@googlegroups.com.

Visit this group at http://groups.google.com/group/methylkit_discussion.
For more options, visit https://groups.google.com/groups/opt_out.

--
You received this message because you are subscribed to the Google Groups "methylkit_discussion" group.

To unsubscribe from this group and stop receiving emails from it, send an email to methylkit_discus...@googlegroups.com.
To post to this group, send email to methylkit_...@googlegroups.com.

[Altuna Akalin]

Hi Kalyan,

I think easiest for you will be option 3. Where you use unite() with min.per.group=10L, then filter the returning object to make sure test group has at least 15 samples covered.

assume obj is the methylBase object. Below is an example of how to achieve this, I haven't tested the code but something along these lines will work.

data=getData(obj)

# this piece gets a TRUE/FALSE data.frame for the test samples
# where TRUE indicates that sample has coverage

filter=!is.na( data[  o...@coverage.index[obj@treatment==0]  ])

# filter based on the number of samples in the test group that have
# at least 15 samples covered 

new.obj=obj[rowSums(filter)>=15,]

but for reference unite function is here:

https://github.com/al2na/methylKit/blob/development/R/backbone.R#L535

[Kalyan K Pasumarthy]

Dear Altuna,

Thanks for the code!

It worked out beautifully. I have tallied it with my version of filter scripted in PERL.

I also learnt some concepts of R within your code!

Regards,
Kalyan

--

[Altuna Akalin]

great! Hope it worked correctly, as I said I didn't test the code.

Best,

Altuna