Excel與基因名的故事（續(xù)）
8號的時候狈癞，刷知乎遇到一個問題履羞。如何評價科學(xué)家重命名了多個人類基因峦萎，以避免被 Excel 自動糾正屡久？沒有別的方法嗎？
雖然不知道知乎從何種途徑觀察到我在研究這個問題然后推了相關(guān)問題給我爱榔。但是這個問題確實推薦的不錯被环。
從問題以及相關(guān)搜索中可以知道，7號详幽，HGNC筛欢，也就是人類基因組組織基因命名委員會，發(fā)了一篇NG的Comment唇聘，介紹了他們新的基因命名指南版姑。
Guidelines for human gene nomenclature 一文中提到：
Symbols that affect data handling and retrieval.
For example, all symbols that autoconverted to dates in Microsoft Excel have been changed (for example, SEPT1 is now SEPTIN1; MARCH1 is now MARCHF1); tRNA synthetase symbols that were also common words have been changed (for example, WARS is now WARS1; CARS is now CARS1).
也就是說，HGNC決定改名了迟郎！
啊剥险，真是喜大普奔。雖然微軟在該事中被詬病宪肖，但是表制，生物學(xué)家還是屈服了。怪不得我之前搜索的時候顯示Sept2的官方名為SEPTIN2匈庭，原來他已經(jīng)改完名了夫凸。
那為什么之前轉(zhuǎn)換成功的是Sept2呢，一方面是數(shù)據(jù)庫沒有及時更新阱持，還有就是，之前是鼠的基因魔熏，而HGNC管理的主要是人的基因衷咽，人的基因名中除了orf之外的字母都是全大寫的。
本篇文章均為人的基因名蒜绽。
續(xù)集
那為什么還有續(xù)集呢镶骗？
因為，我想測試一下代碼的兼容性躲雅，替換sep和mar夠不夠鼎姊，還需不需要替換其他的月份。那么相赁，我就需要一個被Excel影響的基因名清單相寇。
找遍全網(wǎng)居然沒有？
哦不對钮科，有一個唤衫。在公司微信群里，有人提到了一個R包HGNChelper绵脯。它的核心功能就是將不合法的基因名轉(zhuǎn)變?yōu)楹戏ǖ幕蛎ㄒ簿褪荋GNC官方認(rèn)可的）佳励。
稍微介紹一下：
library(HGNChelper)
human = c("FN1", "tp53", "UNKNOWNGENE","7-Sep", "9/7", "1-Mar", "Oct4", "4-Oct",
"OCT4-PG4", "C19ORF71", "C19orf71")
checkGeneSymbols(human)
out:
Human gene symbols should be all upper-case except for the 'orf' in open reading frames. The case of some letters was corrected.x contains non-approved gene symbols x Approved Suggested.Symbol
1 FN1 TRUE FN1
2 tp53 FALSE TP53
3 UNKNOWNGENE FALSE <NA>
4 7-Sep FALSE SEPT7
5 9/7 FALSE SEPT7
6 1-Mar FALSE MARC1 /// MARCH1
7 Oct4 FALSE POU5F1
8 4-Oct FALSE POU5F1
9 OCT4-PG4 FALSE POU5F1P4
10 C19ORF71 FALSE C19orf71
11 C19orf71 TRUE C19orf71
本例來源于：拯救那些被EXCEL篡改的基因名
看起來很厲害休里。但是，考慮到python中Pandas讀取錯誤的基因名讀進(jìn)來后是日期類型赃承，而R讀進(jìn)來是數(shù)字類型妙黍。雖然也不排除其他方式可能能夠讀取Excel中原始的字符2-Sep，但是對我而言瞧剖，這個包并沒有什么用拭嫁。我暫時也不知道他們開發(fā)這個包是用于處理何種來源的原始數(shù)據(jù)。
不過筒繁，這個包的轉(zhuǎn)換方式采用的是存儲數(shù)據(jù)直接映射噩凹，那么，是不是就可以從中提取可能被影響的基因列表呢毡咏。
從他的GitHub下載了rda文件驮宴。
hgnc = new.env()
load('hgnc.table.rda',envir = hgnc)
print(unique(hgnc$hgnc.table$Approved.Symbol)[1:10]) # 數(shù)據(jù)太多，取前10
out:
[1] "NFYC-AS1" "IFITM2" "IFITM3" "PRDX6" "SEC24B-AS1"
[6] "ALDH1L1" "KNOP1" "CYB561D2" "KIR2DL4" "ARHGAP9"
NFYC-AS1也會被影響的嗎？
Alias symbols是0808y08y馆衔，但是這個在Excel中不會被改變啊蚜退。不是很懂。
放棄此路迎罗！
自己動手，豐衣足食
找不到清單片仿，那我就搞一個清單纹安。
基因名，就從最官方的HGNC爬砂豌。
本來打算直接上爬蟲的厢岂，但是仔細(xì)一看，HGNC貼心的提供了HGNC REST web-service阳距。
那就容易了八！！
使用到的主要是兩個接口筐摘。
| 示例 | 類型 | 功能 |
|---|---|---|
| http://rest.genenames.org/fetch/symbol/ZNF3 | fetch | 返回基因名對應(yīng)的信息 |
| http://rest.genenames.org/search/symbol/ZNF3 | search | 搜索基因名卒茬，返回基因列表 |
第二個接口接受類UNIX的匹配符，如ZNF?咖熟、ZNF*圃酵，接口不區(qū)分大小寫。
所有易導(dǎo)致錯誤的基因名都應(yīng)當(dāng)是以月份的簡稱開頭的球恤，所以辜昵，搜索sep2*就可以獲得sep開頭的基因列表，然后拿到每個基因的Previous symbols和Alias symbols咽斧，就可以獲得基因所有的名字了堪置。
首先編寫一個兼容兩個接口的查詢函數(shù)躬存。
import requests
from urllib.parse import urljoin # 用于拼接url
from collections import namedtuple
request_result = namedtuple('request_result',['request_type', 'content']) # 命名元組可以通過屬性訪問值
class RequestTypeError(Exception): # 自定義異常
pass
def request_gene(request_content,request_type = 'search',request_content_type = 'symbol'):
if request_type.lower().startswith('f'):
request_type = 'fetch'
elif request_type.lower().startswith('s'):
request_type = 'search'
else:
raise RequestTypeError # 當(dāng)request_type得小寫不是以f或s開頭時就拋出異常
headers = {'Accept': 'application/json'}
domin = 'http://rest.genenames.org/'
try:
url = urljoin(domin,"/".join([request_type,request_content_type,request_content]))
res = requests.get(url,headers = headers)
content = res.json()
except Exception: # 兼容可能出現(xiàn)的所有錯誤
return None
num_found = content['response']['numFound'] #該內(nèi)容給出結(jié)果的數(shù)量
if num_found:
return request_result(request_type,content['response']['docs']) #docs內(nèi)存儲了具體的內(nèi)容，為列表
else:
return None
然后舀锨，寫循環(huán)執(zhí)行操作岭洲。
def main():
months = ['Jan','Feb','Mar','Apr','May','Jun','Jul','Aug','Sep','Oct','Nov','Dec']
results = []
for month in months:
month_result = request_gene(month+'*',request_type='search') #使用*匹配月份名開頭的基因名
if month_result is None:
continue
for i in month_result.content: #content有多個結(jié)果
gene_symbol = i['symbol']
gene_result = request_gene(gene_symbol,request_type='fetch')
if gene_result is None:
continue
gene_result = gene_result.content[0] # fetch只有一個結(jié)果
if 'alias_symbol' in gene_result:
results.extend([gene_symbol,'alias_symbol',alias_symbol] for alias_symbol in gene_result['alias_symbol'])
if 'prev_symbol' in gene_result:
results.extend([gene_symbol,'prev_symbol',prev_symbol] for prev_symbol in gene_result['prev_symbol'])
return results
獲得結(jié)果，使用Pandas生成Excel坎匿。
results = main()
import pandas as pd
genes = pd.DataFrame(results,columns=['gene_symbol','other_symbol_type','other_symbol'])
genes.to_excel('genes.xlsx',index=False)
打開Excel盾剩。
？替蔬？告私？
說好的自動轉(zhuǎn)換呢〕星牛看來Pandas存儲的Excel在打開時還是能保存原樣的驻粟。不過雙擊單元格再Enter就變過去了，這Excel有點意思凶异。
那怎么辦呢蜀撑，生成csv再打開應(yīng)該就可以了。
genes.to_csv('genes.csv',index=False)
打開剩彬。
果然酷麦，可以了。然后另存為genes_after_open.xlsx
after_open= pd.read_excel('genes_after_open.xlsx')
after_open.head()
整合喉恋，然后篩選出所有被改變的基因名沃饶。
from datetime import datetime
genes['gene_symbol_after_open'] = after_open['gene_symbol']
genes['other_symbol_after_open'] = after_open['other_symbol']
genes[genes.gene_symbol_after_open.apply(lambda x:isinstance(x,datetime))]
官方名字沒有被錯誤更改的（看來改的沒有遺漏）。
再看下other_symbol轻黑。
genes[genes.other_symbol_after_open.apply(lambda x:isinstance(x,datetime))]
out:
看了下都是SEP基因和MARCH開頭的基因绍坝。咦，不是說還有其他月份的嗎苔悦？而且，不是說有一個MARC1和MARCH1被共同轉(zhuǎn)換成1-Mar了嗎椎咧？
查了下玖详，原來MARC1基因的更改后名字居然是MTR開頭的。
那勤讽，看來只進(jìn)行基因名的匹配搜索已經(jīng)不夠了啊蟋座。
看了下接口，官方也提供了全文搜索的選項脚牍，接口類似http://rest.genenames.org/search/ZNF3向臀。
改一下request_gene函數(shù)和main中的調(diào)用方式。
import requests
from urllib.parse import urljoin
from collections import namedtuple
request_result = namedtuple('request_result',['request_type', 'content'])
class RequestTypeError(Exception):
pass
def request_gene(request_content,request_type = 'search',request_content_type = None):
if request_type.startswith('f'):
request_type = 'fetch'
elif request_type.startswith('s'):
request_type = 'search'
else:
raise RequestTypeError
headers = {'Accept': 'application/json'}
domin = 'http://rest.genenames.org/'
try:
if request_content_type is None:
url = urljoin(domin,"/".join([request_type,request_content]))
else:
url = urljoin(domin,"/".join([request_type,request_content_type,request_content]))
res = requests.get(url,headers = headers)
content = res.json()
except Exception:
return None
num_found = content['response']['numFound']
if num_found:
return request_result(request_type,content['response']['docs'])
else:
return None
def main():
months = ['Jan','Feb','Mar','Apr','May','Jun','Jul','Aug','Sep','Oct','Nov','Dec']
results = []
for month in months:
month_result = request_gene(month+'*',request_type='search')
if month_result is None:
continue
for i in month_result.content:
gene_symbol = i['symbol']
gene_result = request_gene(gene_symbol,request_type='fetch',request_content_type='symbol')
if gene_result is None:
continue
gene_result = gene_result.content[0]
if 'alias_symbol' in gene_result:
results.extend([gene_symbol,'alias_symbol',alias_symbol] for alias_symbol in gene_result['alias_symbol'])
if 'prev_symbol' in gene_result:
results.extend([gene_symbol,'prev_symbol',prev_symbol] for prev_symbol in gene_result['prev_symbol'])
return results
results = main()
genes = pd.DataFrame(results,columns=['gene_symbol','other_symbol_type','other_symbol'])
genes.to_csv('genes.csv',index=False)
重新跑诸狭，重復(fù)之前的操作保存生成Excel券膀。
after_open= pd.read_excel('genes_after_open.xlsx')
genes['gene_symbol_after_open'] = after_open['gene_symbol']
genes['other_symbol_after_open'] = after_open['other_symbol']
from datetime import datetime
pd.set_option('display.max_rows', None) #顯示所有行
genes[genes.other_symbol_after_open.apply(lambda x:isinstance(x,datetime))]
看了下君纫，確實包含了之前沒找到的MTAR1基因。
統(tǒng)計一下都有哪些開頭芹彬。
import re
incorrect_symbols = genes[genes.other_symbol_after_open.apply(lambda x:isinstance(x,datetime))].other_symbol
incorrect_symbols.apply(lambda x:re.match(r'([A-Za-z]*)\d*',x).group(1)).value_counts()
out:
SEPT 16
MARCH 11
OCT 9
Mar 9
SEP 6
APR 4
DEC 3
MARC 2
FEB 2
Oct 1
NOV 1
看來確實只有sep和mar需要變換蓄髓。Mar開頭的就離譜∈姘铮看一下会喝。
genes[genes.other_symbol_after_open.apply(lambda x:isinstance(x,datetime))][incorrect_symbols.apply(lambda x:x.lower().startswith('mar')).values]
out:
有一些基因的alias_symbol是Mar開頭的。現(xiàn)在通用的基因名基本都是大寫玩郊，這些應(yīng)該不會再使用了吧肢执，不管他們。MARC和MARCH的問題確實解決不了译红。
最后轉(zhuǎn)換一下预茄，添加一個大寫。
genes['other_symbol_after_transform'] = genes.other_symbol_after_open.apply(lambda x: datetime.strftime(x,'%b%#d').replace('Sep','SEPTIN').replace('Mar','MARCHF').upper() if isinstance(x,datetime) else x)
看一下prev_symbol的復(fù)原情況临庇。
genes[(genes.other_symbol_after_open.apply(lambda x:isinstance(x,datetime))) & (genes.other_symbol_type == 'prev_symbol')]
out:
除了可惡的MARC之外都可以反璃。這個表共31行，應(yīng)該包含了所有27個剛被改名的基因假夺。
再看一下alias_symbol淮蜈。
genes[(genes.other_symbol_after_open.apply(lambda x:isinstance(x,datetime))) & (genes.other_symbol_type == 'alias_symbol')]
out:
這個轉(zhuǎn)換就比較不行了，主要是因為各種牛鬼蛇神都有已卷，SEPTIN6梧田、SEPTIN8居然還有SEP2的alias_symbol，只能在一些symbol上試圖性的恢復(fù)一下侧蘸〔妹校看來有一些基因的alias_symbol也不可靠。
總結(jié)
本文主要是通過找到所有被影響的基因名來驗證之前的替換是否足夠完備讳癌。從結(jié)果來看穿稳，之前的替換已經(jīng)基本夠用，就是需要在更換物種的時候進(jìn)行相應(yīng)的大小寫切換和替換更改晌坤。
