精华内容
下载资源
问答
  • 利用拷贝数畸变,DNA甲基化,基因表达,转录因子表达和miRNA表达数据集,miRDriver应用了基于LASSO的方法来选择miRNA与靶标基因的相互作用。 我们对来自癌症基因组图谱(TCGA)数据库的乳腺癌和卵巢癌数据进行了...
  • 目录Vst介绍计算每个窗口的绝对拷贝数1.文件准备2.编写脚本计算每个窗口的Vst Vst介绍 Vst是通过计算拷贝数的方差来衡量不同群体之间CNV的分化的一个指标,类似于Fst的概念,可以用来鉴定一些高分化的区域。 计算...

    Vst介绍

    Vst是通过计算拷贝数的方差来衡量不同群体之间CNV的分化的一个指标,类似于Fst的概念,可以用来鉴定一些高分化的区域。

    计算方法如下:
    来源:https://doi.org/10.1371/journal.pone.0201611
    Vpop1是指群体1的copy number的方差;
    Vpop2是指群体2的copy number的方差;
    Vtotal是全部个体的copy number的方差;
    Npop指的对应群体的个体数。

    部分文献中,会按照10K的窗口,2K的步长,去计算全基因组上每个窗口的Vst。而一般我们可以用相应的CNV calling的软件检测出每个CNV的位置和拷贝数。因此,我们把这个过程分成两部分:

    1. python脚本计算每个窗口上的绝对拷贝数。
    2. R脚本计算每个窗口的Vst。

    计算每个窗口的绝对拷贝数

    因为每个检测CNV软件生成结果的格式不大相同,这里把结果中CNV的位置和绝对拷贝数整理成一个TSV格式的文件。

    1.文件准备

    在这里插入图片描述

    2.编写脚本

    '''
    @Author: Guo
    @Date: 2020-05-10 10:07:54
    @E-mail: willgyw@126.com
    @Description:  
    @LastEditors: 
    @LastEditTime: 2020-05-13 14:59:18
    '''
    import sys,os
    import logging
    import click
    
    logging.basicConfig(filename=os.path.basename(__file__).replace('.py','.log'),
                        format='%(asctime)s: %(name)s: %(levelname)s: %(message)s',level=logging.DEBUG,filemode='w')
    logging.info(f"The command line is:\n\tpython3 {' '.join(sys.argv)}")
    
    
    @click.command()
    @click.option('-l', '--lengthfile', help='Input length of each chromosome file', required=True)
    @click.option('-t', '--tsv', help='Input copy number tsvfile', required=True)
    @click.option('-w', '--windowsize', type=click.INT, help='windowsize', required=True)
    @click.option('-s', '--stepsize', type=click.INT, help='stepsize', required=True)
    @click.option('-n', '--num', type=click.INT, help='number of samples', required=True)
    @click.option('-o', '--output', help='outputfile',default='copy_number_by_window.tsv')
    
    def main(lengthfile, tsv, windowsize, stepsize, num, output):
        
        step = int(stepsize)
        win = int(windowsize)
        samples_num = int(num)
        fo = open(output, 'w')
        length_dict = {}
        with open(lengthfile, 'r') as lengthf:
            for line in lengthf:
                lines = line.strip().split()
                length_dict[lines[0]] = int(lines[1])
    
        with open(tsv, 'r') as vcf:
            for chr in length_dict:
                for start in range(0, length_dict.get(chr, 0), step):
                    cn_lst = []
                    for num in range(samples_num):
                        cn_lst.append(0)
                    vcf.seek(0)
                    for line in vcf:
                        lines = line.strip().split()
                        c = lines[0]
                        s = lines[1]
                        e = lines[2]
                        if c == chr:
                            if int(start) < int(e) < (int(start) + int(win)):
                                for num in range(samples_num):
                                    cn_lst[num] += round(float(lines[(num+3)].split(':')[1]), 2)
                            else:
                                continue
                        else:
                            continue
                    if sum(cn_lst) > 0:
                        cn_lst_new = [str(round(x, 2)) for x in cn_lst]
                        fo.write('\t'.join([str(chr), str(start), str(int(start) + int(win))]) + '\t' + '\t'.join(cn_lst_new) + '\n')
        fo.close()
    
    if __name__ == '__main__':
        main()
    
    

    脚本利用click模块引入命令行参数,需要提供每个染色体的长度、整理好的tsv格式文件、窗口、步长、样本数等参数。

    最后会生成类似下图的文件:

    结果文件
    窗口位置和每个窗口的拷贝数。

    计算每个窗口的Vst

    有了上述的结果文件之后,我们就可以利用R脚本计算每个窗口的Vst了。

    R的脚本我就不放了,因为不是我写的,尊重原创。python脚本还有很多可以完善的地方,欢迎大家批评指正。

    展开全文
  • CoDeCZ从基因或外显子水平的靶向测序... 结果不是真正的拷贝数,但是是半定量的:它使用修改后的z分数来计算每个指定区域的标准化覆盖率的偏差。 Z分数越高,副本数越高。 最终的z分数取决于数据的质量和参考池的组成。
  • 拷贝数变化模拟器(CNV-Sim) 在基因组学中,拷贝数变异(CNV)是基因组中一种结构变异,其中基因组的各个部分被复制或缺失。 人群中个体之间变异(重复/缺失)的数量有所不同。 拷贝数变异模拟器(CNV-Sim)是一...
  • 没使用ccur之前拷贝过来的是科学计数法的字符串。使用后是字符串只能计数不能计算。继续解决 "Update [项目跟踪表$],[backlog$] set [项目跟踪表$].backlog = ccur([backlog$].backlog) where [项目跟踪表$].psaid=...

    "Update [项目跟踪表$],[backlog$] set [项目跟踪表$].backlog = ccur([backlog$].backlog) where [项目跟踪表$].psaid=[backlog$].项目编号"

    没使用ccur之前拷贝过来的是科学计数法的字符串。使用后是字符串只能计数不能计算。继续解决


    alter  table  [项目跟踪表$]   alter column backlog float

    以上语句无效,后来查了下excel没有数据类型限制,执行上面的报错。


    执行几次后出现:

    数据库引擎已停止进程,因为您和其他用户试图同时改变同一数据 

    想想我们用JAVA的时候都用try catch finally在异常是关闭连接。VBA里没这个东东。可能是连接没关闭导致。


    补充异常处理:

     On Error GoTo Err_Handle

      ……

    Err_Handle:
        MsgBox Err.Description
        cnn.Close
        Set cnn = Nothing


    虽然加了这个东东,前面的引擎停止仍提示中,网上各种解决办法还要下什么工具。就是写个小工具弄这个太烦了。

    终于在微软官网看到把可能引擎冲突的数据删掉就可以了。果断拷贝删除。好了。


    继续解决拷贝过来的是字符串不能计算的问题。

    展开全文
  • 全基因组范围内的绝对拷贝数和拷贝中性变异的全面可视化 联系人: Victor Renault / Alexandre How-Kit() aCNV​​iewer(Absolute CNV Viewer)是一种工具,可以在一组癌症样本中可视化绝对CNV和cn-LOH。 aCNV​...
  • Python的优势目前应用特点Part 2 6种内置对象1 整数、浮点数2 字符串3 列表4 元组5 字典6 集合什么是浅拷贝shallow copy案例2.1 编写程序,根据输入的半径,计算圆的面积案例2.2 编写程序,实现凯撒密码,按照规定的...

    为什么写这篇文章

    这是第一天学习python ,由于本科基础不太好,又是学习新知识内容、故决定从写文章开始做到一定量的输出,以便对知识进行巩固、消化吸收。同时也希望给希望快速学习了解python的朋友一点帮助、小弟水平有限,文中有理解不对,不到位的地方,还请指正!

    Part 1 入门部分

    什么是编程语言?

    编程语言也是一门语言,它在面向人类所熟知的自然语言的语义空间同面向机器的机器语义空间架起了一座桥梁,进行映射、
    从最早的机器语言Machine Language,即二进制代码发展到了汇编语言Assembly Language,封装了诸多关键字,提升了一定可读性,但汇编语言与机器语言仍是一一对应的。虽效率高,但若进行大规模编程,那将是极为痛苦的、故发展出了高级语言High level Language,如图C/C++,JAVA,PYHON等流行的语言。高级语言不论从各维度都更符合人类学习,符合建模坚决问题的思维习惯,满足应用场景的需求、

    Python的优势

    选择Python是因为其用途相对广泛、且简单易学

    目前应用

    从官网中:

    Python’s convenience has made it the most popular language for machine learning and artificial intelligence. Python’s flexibility has allowed Anyscale to make ML/AI scalable from laptops to clusters.(https://www.python.org/)

    可见Python对人工智能领域和机器学习领域是十分便捷的、能结合本专业使用到的方向有网络爬虫、机器学习、数据分析、神经网络等

    特点

    开源、简单、不断进化

    Part 2 6种内置对象

    1 整数、浮点数

    整数int顾名思义,0,1,2,3这些没有小数点的都是整数范畴、
    浮点数float就是带有小数点的数,之所以叫浮点数,是因为小数点可通过科学计数法表示时浮动、

    python不需要像C++,java制定类型声明变量,如int a=0;int b=0.3;
    可以直接赋值,如a=3 a=3.0
    结尾也不需要; 封号
    但是如果要浮点数的3.0也可以声明float(3),输出即为3.0,十分方便

    需要掌握四则运算+ -;* / 和 % ;divmod ()得到除数和被除数;round() 可以四舍五入保留小数; type()可以获取类型;
    发现计算中存在的诸多异常,并通过标准库中的math和其他工具解决之
    如以下异常

    // 异常
    >>> 0.1+0.2
    0.30000000000000004//问题输出
    >>> round(0.1+0.2,1)
    0.3//调整后输出
    

    可以这样解决:(当然也有很多别的解决方法)

    >>> import decimal
    >>> a=decimal.Decimal('0.1')
    >>> b=decimal.Decimal('0.2')
    >>> a+b
    Decimal('0.3')
    

    需要掌握
    dir([object]) 可以返回模块的属性列表(查看对象的属性和方法)
    help([object])的用法 可以返回函数和模块用途的详情说明!

    2 字符串

    字母、标点符号、控制符都是字符
    字符编码: ACSII、Unicode、UTF-8
    使用引号(“”或’’)创建字符串即可
    需要在字符串中使用特殊字符时可以用反斜杠(\)转义字符,如下表示:

    >>> b='This is YM. I took my brother's bag'//识别错误前两个''为一对
      File "<stdin>", line 1
        b='This is YM. I took my brother's bag'
                                         ^
    SyntaxError: invalid syntax
    >>> b='This is YM. I took my brother\'s bag'
    >>> b
    "This is YM. I took my brother's bag"//成功
    

    索引(左右开始)

    >>> help(str.index)
    Help on method_descriptor:
    index(...)
        S.index(sub[, start[, end]]) -> int
    
        Return the lowest index in S where substring sub is found,
        such that sub is contained within S[start:end].  Optional
        arguments start and end are interpreted as in slice notation.
    
        Raises ValueError when the substring is not found.
        >>> s="I love python"
    >>> s
    'I love python'
    >>> s.index('love')
    2
    >>> s.index('o',6)
    11
    

    切片

    >>> s.split(" ")//拆分成列表
    ['I', 'love', 'python']
    //
    >>> q="good"
    >>> q[1:3]//按照下表为12切片
    'oo'
    python
    //
    >>> q
    'abcdefg'
    >>> q[2:5]
    'cde'//切片2-4位(第一位是0

    切片再通过制定字符连接

    >>> lst=s.split(" ")
    >>> lst
    ['I', 'love', 'python']
    >>> "~".join(lst)
    'I~love~python'
    

    此外,学习了两个内置函数,input() output()、外加字符串格式化输出format()

    3 列表

    顺口溜:列表是个筐,什么都能装。
    列表是一个序列,也是容器、其中的元素是可变的。
    与字符串相同地存在索引(也可从右往左,即-1开始)

    >>> p
    ['a', 'b', 1, 'd']
    >>> p[0]
    'a'
    >>> p[0]=2
    >>> p
    [2, 'b', 1, 'd']
    

    和 切片,

    >>> p
    ['a', 'b', 1, 'd']
    >>> p[1:3]
    ['b', 1]
    

    方法差不多在此不赘述

    列表和字符串不同的地方:不可以按照索引修改值

    >>> q="abcdefg"
    >>> q[2]=w
    Traceback (most recent call last):
      File "<stdin>", line 1, in <module>
    NameError: name 'w' is not defined
       //字符串按照索引改值失败
    >>> p=['a','b','c','d']
    >>> p[2]=1
    >>> p
    ['a', 'b', 1, 'd']
    //相反,列表成功
    

    其余均相同,如加法,乘法,len() ,in 等

    列表中增加元素有append(),insert(),extend(),注意何为可迭代对象、多参照帮助文档
    删除元素有 remove(),pop(),clear()
    排序sort() 反序reverse() 等等、需要用的时候在了解

    与字符串相互转化list(x); “”.join(lst)

    >>> st='this is a test string'
    >>> st
    'this is a test string'
    >>> lst=list(st)
    >>> lst//成功从string转到list
    ['t', 'h', 'i', 's', ' ', 'i', 's', ' ', 'a', ' ', 't', 'e', 's', 't', ' ', 's', 't', 'r', 'i', 'n', 'g']
    >>> "".join(lst)
    'this is a test string'//list转到string
    

    4 元组

    需要注意的是,只有一个元素的元组这样定义tuple(‘a’,)

    >>> b=tuple('a',)
    >>> b
    ('a',)
    >>> c=['a']
    >>> c
    ['a']
    

    元组是一个序列,容器、区别是在于其中的元素不可修改值,也带来了其优点:运算速度快、

    索引(正向,反向)和切片

    基本操作+、*、len()、in与列表基本相同,不多做介绍

    元组与列表相互转换lst=list(X)/u=tuple(X)

    >>> t=('ts',1,2,[1,2,3])
    >>> t
    ('ts', 1, 2, [1, 2, 3])
    >>> lst=list(t)
    >>> lst
    ['ts', 1, 2, [1, 2, 3]]
    >>> u=tuple(lst)
    >>> u
    ('ts', 1, 2, [1, 2, 3])
    >>>
    

    5 字典

    字典是一组组映射关系,也可以理解为键值对 ,即 key-value(key不可重复)

    d = {key1 : value1, key2 : value2 }

    读d.get(‘a’,‘lanqiao’)

    其他/setdefault/update/pop/popitem

    字典和列表:字典不是序列,列表是序列

    两者都是容器类对象,可变对象,从Python3.6开始,字典也有顺序

    一个例题

    写函数,检查传入字典的每一个value的长度,如果大于2那么仅保留前两个长度的内容,并将新内容返回给调用者:PS字典中的value只能是字符串或列表

    def check(**kwargs):
        n={}
        for k,v in kwargs.items():
             if len(v)>2:#筛子
                  v=v[:2]
             n[k]=v
        return n
    
    empty={}
    flag=1
    while(flag):
         i=input("k-v:")
         if i=='exit':
             flag=0
             print("exit")
         else:
             a,b=map(str,i.split('-'))
             empty[a]=b
             print(empty)
             
    s=check(**empty)
    print(s)
    

    测试:

    k-v:tom-basketball
    {'tom': 'basketball'}
    k-v:alice-23234
    {'tom': 'basketball', 'alice': '23234'}
    k-v:bob-23jiojf
    {'tom': 'basketball', 'alice': '23234', 'bob': '23jiojf'}
    k-v:exit
    exit
    {'tom': 'ba', 'alice': '23', 'bob': '23'}
    

    6 集合

    集合是一个拥有确定(唯一)的,且元素无序的可变的数据组织形式

    分为可变集合(集合中的元素可以动态的增加或删除–set)和非可变集合(集合中的元素不可改变–frozenset)

    可变集合的定义:

    >>> s=set('aabbcded')
    >>> s
    {'d', 'c', 'b', 'a', 'e'}
    

    add,pop,remove,discard

    不可变集合的定义:

    >>> ss=frozenset('qiwsir')
    >>> ss
    frozenset({'r', 's', 'q', 'i', 'w'})
    

    集合的特点:里头的对象应该是不可变对象、列表中是可变的。即集合中不能包含列表

    关系和运算:is superser() is superset(b) sisubset(a) union 交并补差集等等

    什么是浅拷贝shallow copy

    >>> help(list.copy)
    Help on method_descriptor:
    
    copy(self, /)
        Return a shallow copy of the list.
    
    >>> help(dict.copy)
    Help on method_descriptor:
    
    copy(...)
        D.copy() -> a shallow copy of D
    
    >>> help(set.copy)
    Help on method_descriptor:
    
    copy(...)
        Return a shallow copy of a set.
    
    >>>
    

    可以发现,观察帮助文档,三种容器(列表、字典、集合)所提供的copy()方法都是做的shallow copy,也就是浅拷贝,那么到底什么是浅拷贝?先给出结论:浅在哪里?拷贝的是第一层、里面的层并未拷贝
    以列表为例,观察如下代码:

    >>> lst1=['brian',19,['java','python','C++']]
    >>> lst1
    ['brian', 19, ['java', 'python', 'C++']]
    >>> type(lst1)
    <class 'list'>
    >>> lst2=lst1.copy()
    >>> lst2
    ['brian', 19, ['java', 'python', 'C++']]
    

    首先创建了一个lst1,并通过type()查看到了其中类型为list
    通过列表中的拷贝函数、拷贝出了lst2,可以发现显示结果一样

    >>> id(lst1)
    1633969095872
    >>> id(lst2)
    1633964907200
    >>> lst1 is lst2
    False
    >>> lst1[0] is lst2[0]
    True
    >>> lst1[2] is lst2[2]
    True
    

    观察lst1和lst2的地址,是不同的
    且通过 is 来比较,得知lst1和lst不是一个对象

    >>> lst1[0]=100
    >>> lst1
    [100, 19, ['java', 'python', 'C++']]
    >>> lst2
    ['brian', 19, ['java', 'python', 'C++']]
    

    这个时候对第一层的内容进行修改,修改lst1中的第一个元素为100,发现lst2中第一个元素仍然是brian
    可以知第一层确实不存在引用关系

    >>> lst1[2][0]='C#'
    >>> lst1
    [100, 19, ['C#', 'python', 'C++']]
    >>> lst2
    ['brian', 19, ['C#', 'python', 'C++']]
    

    当把lst1中的第一个元素java改为C#后
    再输出lst2,发现lst2中的[2][0]元素也被改变了、故这就是浅拷贝
    如何解决呢、可以通过 import copy、调用deepcopy()来实现深拷贝

    >>> import copy
    >>> lst3=copy.deepcopy(lst1)
    >>> lst3
    [100, 19, ['C#', 'python', 'C++']]
    >>> lst1[2][0]='VB'
    >>> lst1
    [100, 19, ['VB', 'python', 'C++']]
    >>> lst3
    [100, 19, ['C#', 'python', 'C++']]
    

    如上代码,这时你就会发现,lst3与lst1不存在引用了,即消除了浅拷贝的隐患

    案例2.1 编写程序,根据输入的半径,计算圆的面积

    import math
    radius=(float)(input("请输入半径"))
    res=math.pi*radius*radius
    print(round(res,3))
    //运行结果
    请输入半径2
    12.566
    

    案例2.2 编写程序,实现凯撒密码,按照规定的偏移位加密

    利用凯撒密码方案,实现对用户输入文字的加密操作
    凯撒密码(Caesar cipher),是一种替换加密的技术,明文中的所有字母都在字母表上向后(向前)
    按照一个固定数目进行偏移后背替换成密文。如:偏移量是3时,所有的子母A被替换成D,B变成E
    由于知识结构缺陷,先做单个字母
    s=str(input("input a letter:"))
    n=3
    o=chr(ord(s)+3)
    print(s,"后移动",n,"位之后是:",o)
    //测试结果
    input a letter:a
    a 后移动 3 位之后是: d
    

    案例2.3 编写程序,实现英文大小写转换

    s=input("输入英文:")
    lst=[]
    for i in s:
        if i.islower():
            lst.append(i.upper())
        else:
            lst.append(i.lower())
    n="".join(lst)
    print(s,"==>",n)
    
    //结果
    输入英文:Apple
    Apple ==> aPPLE
    

    作业一 用户输入国家名称,打印出所输入的国家名称和首都

    a=input("input country:")
    d={'china':'beijing','japan':'tokyo','russia':'莫斯科','england':'london','german':'柏林'}
    print(a,"\'s captial is:",d.get(a))
    //结果
    input country:china
    china 's captial is: beijing
    

    作业二 用户输入数字,显示对应的英文

    a=input("input number:")
    lst=list(a)
    re=[]
    d={'1':'one','2':'two','3':'three','4':'four','5':'five','6':'six','7':'seven','8':'eight','9':'nine','0':'zero'}
    for i in lst:
        re.append(d.get(i))
    n=" ".join(re)
    print(n)
    //结果
    input number:1024
    one zero two four
    

    Part 3 基本编程语句

    布尔类型

    True 1 真
    False 0 假

    比较运算

    比较运算就是比较两个对象的大小
    如 >、 <、 ==、 !=、 >=、 <=

    逻辑运算(布尔运算)

    逻辑运算多用于复合条件的判断 and or not

    >>> 3>2 and 3>1
    True
    >>> 3<2 and 3>10
    False
    >>> 1 or 0
    1
    >>> 0 or 0
    0
    >>> 1 or 1
    1
    >>> not 3
    False
    >>> not 0
    True
    

    赋值语句

    基本形式:variable =object
    其他花样:

    >>> a=1,2,3
    >>> a
    (1, 2, 3)
    //first
    >>> a,b,c=1,2,3
    >>> a
    1
    >>> b
    2
    >>> c
    3
    //second
    >>> a,_,c=1,2,3
    >>> _
    2
    >>> c
    3
    //third
    >>> a,b,*c=1,2,3,4,5
    >>> a
    1
    >>> b
    2
    >>> c
    [3, 4, 5]
    //forth
    >>> a=b=1
    >>> b
    1
    >>> a
    1
    //fifth
    

    条件语句

    不同之处是不需要括号,通过:和缩进来判断

    if(<expr>):
        <statment>
    elif(<expr>):
        <statement>
    else:
        <statement>
    <following_statment>
    

    for循环语句

    for(循环规则):
    [ 空四格 ]语句块

    常用函数: range() zip() enumerate() 列表解析

    while循环语句

    while [condition]:
    [ 空四格 ]statements

    小例题:创建一个数据及,包含1到10的随机整数,总共100个数。统计每个数字的次数

    import random
    lst=[]
    for i in range(100):
        n=random.randint(1,10)
        lst.append(n)
    
    d={}
    for n in lst:
        if n in d:
            d[n]+=1
        else:
            d[n]=1
    print(d)
    //结果如下
    {2: 9, 4: 6, 3: 14, 7: 8, 9: 12, 1: 7, 8: 12, 5: 13, 10: 13, 6: 6}
    

    案例3.1创建一个列表,其中的元素是100以内的能被3整除的正整数

    r=[]
    for i in range(100):
        if ((i+1)%3)==0:
            r.append(i+1)
    print(r)
    //结果
    [3, 6, 9, 12, 15, 18, 21, 24, 27, 30, 33, 36, 39, 42, 45, 48, 51, 54, 57, 60, 63, 66, 69, 72, 75, 78, 81, 84, 87, 90, 93, 96, 99]
    PS C:\Users\ASUS\Desktop\MLCode\python>
    

    案例3.2 字符串 s=‘Life is short You need python’。统计这个字符串中每个单词的字母数量

    s='Life is short You need python'
    d={}
    for i in s:
        if i.isalpha():
            if i in d:
                d[i]+=1
            else:
                d[i]=1
    print(d)
    //结果
    {'L': 1, 'i': 2, 'f': 1, 'e': 3, 's': 2, 'h': 2, 'o': 3, 'r': 1, 't': 2, 'Y': 1, 'u': 1, 'n': 2, 'd': 1, 'p': 1, 'y': 1}
    

    案例3.3 制作一个满足下功能的猜数游戏

    计算机随机生成一个100以内的正整数;
    用户通过键盘输入数字,猜测计算机所生成的随机数;
    注意:对用户的输入次数不做限制;

    import random
    m=random.randint(1,100)
    flag=1
    while(flag):
        a=(int)(input("please input your num from 1 to 100:"))
        if a==m:
            print("Success!")
            flag=0
        elif(a<m):
            print("smaller than res")
        else:
            print("bigger than res")
    print("gamer over")
    //结果
    please input your num from 1 to 100:50
    bigger than res
    please input your num from 1 to 100:25
    smaller than res
    please input your num from 1 to 100:37
    bigger than res
    please input your num from 1 to 100:30
    bigger than res
    please input your num from 1 to 100:28
    bigger than res
    please input your num from 1 to 100:26
    Success!
    gamer over
    

    一下学的有点多、就先到这里~ 未完待续

    展开全文
  • 今天我们学习一个拷贝数变异的整合软件——GISTIC2。注意,这和软件本身并不做CNV calling,而是主要用于检测一组样品中显着扩增或缺失的基因组区域(明白一点说就是你需要提供一批样本中的每个样本的CNV检测结果,...

    今天我们学习一个拷贝数变异的整合软件——GISTIC2。注意,这和软件本身并不做CNV calling,而是主要用于检测一组样品中显着扩增或缺失的基因组区域(明白一点说就是你需要提供一批样本中的每个样本的CNV检测结果,这个软件经过呼啦呼啦显著性计算会告诉你这一批样本中显著扩增和缺失的是哪些区域)。这个是癌症基因组CNV分析中十分常见也十分必要的内容。

     

    1.软件安装

    注意:

    a.软件包没有打包在一个文件夹下,所以第2步新建了一个GISTIC2文件夹,请在该文件夹下解压源文件;

    b.第5、6步安装Matlab环境的时候,-destinationFolder后面接的是绝对路径,请将“PATH”替换为GISTIC2文件夹所在的绝对路径。

    c.下面步骤编译好后/PATH/GISTIC2/gistic2即为可执行文件。

    wget -c ftp://ftp.broadinstitute.org/pub/GISTIC2.0/GISTIC_2_0_23.tar.gz
    mkdir GISTIC2 && mv GISTIC_2_0_23.tar.gz GISTIC2/ && cd GISTIC2/
    tar -xzvf GISTIC_2_0_23.tar.gz
    cd MCR_Installer/ && unzip MCRInstaller.zip
    ./install -mode silent -agreeToLicense yes -destinationFolder /PATH/GISTIC2/MATLAB_Compiler_Runtime/
    export XAPPLRESDIR=/PATH/GISTIC2/MATLAB_Compiler_Runtime/v83/X11/app-defaults:$XAPPLRESDIR

     

    2.输入文件准备

    a.Segmentation File

    文件应含有所有样本的CNV信息,共6列,以tab键分割:

    (1)  Sample:样本名称

    (2)  Chromosome:染色体名称

    (3)  Start Position:起始位置

    (4)  End Position :终止位置

    (5)  Num markers:在区域内的标记数目

    (6)  Seg.CN:拷贝数取log2()-1

    文件示例:

    b.Reference Genome File

    参考基因组文件,软件包里的refgenefiles文件夹下提供了不同版本的参考基因组mat文件,这里我是用的hg19参考基因组,因此选择hg19.mat文件作为输入的参考基因组文件。

     

    3.软件使用

    软件的参数有很多,详细可以见软件包下的document,这里介绍几个常用的:

    -b  输出目录,后面接一个输出文件前缀prefix

    -seg  上面讲到的Segmentation文件

    -refgene  上面讲到的参考基因组文件

    -conf  置信水平,默认是0.75

    gistic2  -b /outdir/prefix -seg Segmentation_File -refgene hg19.mat -conf 0.9

     

    4.结果文件解读

    运行结束后,会在输出目录下生成一些文件和图片:

    a.all_lesions.conf_XX.txt(XX是置信度)

    该文件是gistic全部结果的整合,包含了扩增和缺失区域,以及每个区域中扩增或缺失来源于哪些样本。从第10列开始往后是单样本信息,下面解释前9列:

    (1)  Unique Name:  鉴定出的扩增或缺失区域名称;

    (2)  Descriptor:  位于基因组染色体臂的位置;

    (3)  Wide Peak Limits:  最可能包含目标基因的参考基因组边界范围;

    (4)  Peak Limits:  最大扩增或缺失区域的参考基因组边界范围;

    (5)  Region Limits:  显著扩增或缺失区域的参考基因组边界范围;

    (6)  q-values:  峰值区域的q值;

    (7)  Residual q-values:  去除与相同染色体中其他更显著的峰区域重叠的扩增或缺失后,峰区域的q值;

    (8) Broad or Focal:  根据置信度将事件分为“broad”和“focal”;

    (9) Amplitude Threshold:  文件10列及以后列的每个样本给出值的含义解释;

    b.amp_genes.conf_XX.txt(XX是置信度)

    该文件是扩增区域的结果,每一列是一个扩增峰的结果,文件共4行:

    (1)  cytoband

    (2)  q-value

    (3)  residual q-value

    (4)  wide peak boundaries

    c.del_genes.conf_XX.txt(XX是置信度)

    该文件是缺失区域的结果,每一列是一个缺失峰的结果,文件共4行,格式和Amplification Genes File相同。

    d.Gistic Scores File

    GISTIC打分结果文件,文件列出了整个基因组中扩增和缺失的q值[表示为-log10(q)],G得分,异常样本之间的平均幅度以及畸变频率。文件可以导入到IGV中。

    e.raw_copy_number.pdf

    一个pdf热图文件, 基因组沿垂直轴表示,样品水平排列,红色为拷贝数扩增,蓝色为拷贝数缺失。

    f.amp_qplot.pdf、del_qplot.pdf

    扩增或缺失区域单独的展示。

    更多生信小知识关注:

    展开全文
  • Whinhin-SamplE拷贝数畸变检测仪 使用全基因组数据检测母体血浆样品中的胎儿三体性和较小的CNV。 作者: 本文描述了所使用算法的: 。 当前版本未在论文中单独发布。 更改列表 此版本的WISECONDOR与以前发布的版本...

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 857
精华内容 342
关键字:

拷贝数计算