下載基因組外顯子信心
網站 ftp://ftp.ncbi.nlm.nih.gov/pub/CCDS/current_human/
wget ftp://ftp.ncbi.nlm.nih.gov/pub/CCDS/current_human/CCDS.current.txt
運行下列代碼 得到外顯子大約36M
import
re
import
os
from
collections
import
OrderedDict
from
operator
import
itemgetter
os.chdir(
'/Users/yangqin/Desktop/biotree/python/'
)
exonLength
=
0
# 外顯子長度尚未計算
overlapExons
=
OrderedDict()
with
open
(
'CCDS.current.txt'
,
'rt'
) as f:
# 把CCCDS.current.txt打開定義為f
for
line
in
f:
# 逐行讀取f
if
line.startswith(
'#'
): # 如果該行以井號開頭,該行為表頭
continue
# 跳過該行
line
=
line.rstrip()
# 把右邊的空格去掉,_.rstrip()把右邊的空格去掉
lst
=
line.split(
'\t'
)
# 用split分開該行的各列數據,以\t分隔
if
lst[
-
2
]
=
=
'-'
:
# 如果倒數第二列是‘-’,沒有外顯子
continue
# 跳過該行
lst[
-
2
]
=
re.sub(
'\[|\]'
,
' '
, lst[
-
2
])
# 正則表達式,去掉倒數第二列兩端的[],換成空格
exons
=
lst[
-
2
].split(
', '
)
# 把該列(每個基因)多個外顯子隔開
for
exon
in
exons:
# 對於每一個外顯子
start
=
int
(exon.split(
'-'
)[
0
])
# 拿到該外顯子的起始坐標
end
=
int
(exon.split(
'-'
)[
1
])
# 拿到該外顯子的終止坐標
coordinate
=
lst[
0
]
+
':'
+
exon
# 定義每一個外顯子,以防重復計算
if
coordinate
not
in
overlapExons.keys():
# 如果該外顯子尚未計算
overlapExons[coordinate]
=
1
# 加入該外顯子
exonLength
+
=
end
-
start
# 計算長度
print
(exonLength)
36061455