Chinaunix

标题: 计算投影总长度 [打印本页]
作者: junlingpang    时间: 2011-06-09 20:37
标题: 计算投影总长度
A    5    10
A    15    18
A    2    7
A    4    6
B    1    10
B    3    5
C    1    3
C    4    5
C    11    13

如上，第一列为标签（列间以tab分隔），其后两列可以看成是x轴上一条线段的两个端点，计算相同标签下线段的投影总长度（即不计算重叠的部分）。

以标签A为例：

最终A标签下线段投影总长度为10-2+18-15=11。

请大侠们指点一下怎么编程序解决吧。

谢谢啦^_^
作者: yybmsrs    时间: 2011-06-09 21:44
分情况一个个算？
作者: junlingpang    时间: 2011-06-09 21:49
回复 2# yybmsrs


    嗯，得考虑到每一种情况，思路有点儿乱，谢谢参与
作者: iamlimeng    时间: 2011-06-09 22:15
没严格测试，你看看是否正确。
  1. #!/usr/bin/perl

  3. use strict;
  4. use warnings;

  6. my %result;
  7. while (<DATA>) {
  8.         my ($key,$start,$end) = split /\s+/;
  9.          foreach ($start+1 .. $end) {
  10.                  $result{$key}{'count'}++ unless ($result{$key}{$_});
  11.                  $result{$key}{$_} = 1;
  12.          }
  13. }
  14. foreach (keys %result) {
  15.         print "$_\t$result{$_}{'count'}\n";
  16. }

  18. <STDIN>;

  20. __DATA__
  21. A    5    10
  22. A    15    18
  23. A    2    7
  24. A    4    6
  25. B    1    10
  26. B    3    5
  27. C    1    3
  28. C    4    5
  29. C    11    13
复制代码
作者: zhlong8    时间: 2011-06-09 22:29
A 并 B = A + B - (A 交 B) 这个公式叫什么来着？可以用
作者: junlingpang    时间: 2011-06-09 22:31
回复 4# iamlimeng


    刚试了下，可以的，谢谢了
    我再好好研究研究，呵呵
作者: junlingpang    时间: 2011-06-09 22:34
回复 5# zhlong8


    嗯。对于两个的情况还行，多个并集可能实现起来有点儿麻烦。谢谢提供思路
作者: zhlong8    时间: 2011-06-09 22:36
回复  zhlong8


    嗯。对于两个的情况还行，多个并集可能实现起来有点儿麻烦。谢谢提供思路
junlingpang 发表于 2011-06-09 22:34



    对计算机来说多个和两个没有质的区别，它就是做重复运算的。而且写出来的程序具有通用性，抽象级别很高
作者: junlingpang    时间: 2011-06-09 22:41
回复 8# zhlong8


    嗯，我试一下，谢谢大侠
作者: zhlong8    时间: 2011-06-09 22:52
回复 9# junlingpang


    这是最简单的递归思路，其实数端点的方法更快
作者: junlingpang    时间: 2011-06-09 22:58
回复 10# zhlong8


    大侠，能具体讲一下数端点的方法吗？我刚学编程，急需指点，谢了
作者: zhlong8    时间: 2011-06-09 23:08
回复  zhlong8


    大侠，能具体讲一下数端点的方法吗？我刚学编程，急需指点，谢了
junlingpang 发表于 2011-06-09 22:58



    将这些点排好序，如果开始端点和结束端点数相等即为一段连续区间（当然是最短优先）所以只需要数一次就可以知道总长度了
作者: zhlong8    时间: 2011-06-09 23:15
本帖最后由 zhlong8 于 2011-06-09 23:17 编辑

用 scheme 写个递归的代码
; 只需要定义两个基本元素的交 inter
; 单个元素的值 value
; 用到两个公式
; A ∪ B = A + B - A ∩ B
; A ∩ (B1 ∪ B2 ... Bn) = (A ∩ B1) ∪ (A ∩ B2) ... ∪ (A ∩ Bn)
(define (unions items)
  (if (null? (cdr iems))
    (value (car items))
    (let ((a (car items))
          (b (cdr items)))
      (+ (value a) (unions b)
         (- (unions (map (lambda (x) (inter a x)) b)))))))
作者: junlingpang    时间: 2011-06-09 23:26
回复 13# zhlong8


    啊，谢谢了，我得好好消化一下，多谢多谢^_^
作者: chenhao392    时间: 2011-06-10 02:27
不知道楼主是想做什么，生物信息里面有merge gene reads的工具(C binary)，你可以假装是ABC三个染色体，这些片段是gene的位置，直接拿到结果，如有兴趣，可回复此帖联系我
作者: guap514    时间: 2011-06-10 08:11
  1. use List::Util qw(sum);
  2. my $hash;
  3. while(<DATA>){
  4.     my @t = split /\s+/;
  5.     $hash->{ $t[0] }->[ $_ ] = 1 for $t[1] + 1 .. $t[2];
  6. }

  8. foreach my $k ( keys %$hash ){
  9.     print "$k => ", sum( @{ $hash->{ $k } } ), "\n";
  10. }

  12. __DATA__
  13. A    5     10
  14. A    15   18
  15. A    2    7
  16. A    4    6
  17. B    1    10
  18. B    3    5
  19. C    1    3
  20. C    4    5
  21. C    11    13
复制代码
作者: junlingpang    时间: 2011-06-10 09:01
回复 15# chenhao392


    啊~太好啦，呵呵，我也是学生物信息的，研一，刚入门~
  最近在做关于基因中重复序列的一个小的实验，拿到的基因组中有标注基因的位置和重复序列的位置信息，我想统计一下基因中重复序列的百分比。因为一个基因中往往不止一个重复序列，就像上面举的例子那样，所以一开始想编个程序解决一下。
  您说的那个软件可以解决类似的问题吗？谢谢指点^_^
作者: junlingpang    时间: 2011-06-10 09:07
回复 16# guap514


    和4楼的思路类似，但是您的代码更简洁，谢谢大侠帮忙^_^
作者: chenhao392    时间: 2011-06-10 09:26
回复 17# junlingpang


    呵呵，确实可以利用这个工具，不可避免的要用到Bed格式，不过既然都做序列了，这个格式也是经常用的吧
    请google Bedtools manual， 找到里面关于mergeBed/coverageBed 的章节阅读。


    如果只想算百分比，就直接用coverageBed，如果还想 继续这个topic，即 merge这些location再算除法，就看mergeBed.

    咱们干这行的，能利用别人工具是效率最高的了，什么Samtools啊，Bedtools啊，都可以瞅瞅..
作者: junlingpang    时间: 2011-06-10 09:32
回复 19# chenhao392


    嗯，太好了，多谢多谢，刚起步要学习的还有好多哇~谢喽
作者: 2010chris    时间: 2011-06-13 12:36
学习了



欢迎光临 Chinaunix (http://bbs.chinaunix.net/) Powered by Discuz! X3.2