遺伝子の起源 -

生物が持つ遺伝子は、翻訳されて、数百〜数千アミノ酸のタンパク質になるのが一般的である。数百〜数千残基の長大なタンパク質or遺伝子が、ランダムにアミノ酸やヌクレオチドの重合を繰り返しただけで出現する確率は、普通に考えると極めて低い。ヒトの遺伝子で最も巨大なものの一つはtitin遺伝子で、363個のエクソンで構成され、最大のisoformは36000残基近くになるようである。

一般的に、生物の持つタンパク質は、複数の共通するドメインからなることが知られていて、ヒトの遺伝子数は2〜3万とされている(選択的スプライシングがあるので翻訳されるタンパク質の種類は10万種類くらいになると考えられている)が、1000個程度のドメインの組み合わせでバリエーションを増やしているらしい。また、真核生物、真正細菌、古細菌間でも、共通するドメインの利用が確認されているとされる。

cf)Evidence for PDZ domains in bacteria, yeast, and plants.
https://dx.doi.org/10.1002%2Fpro.5560060225

SMART: Domain List
http://smart.embl-heidelberg.de/smart/domain_table.cgi
を見ると、1400以上のドメインが列挙されている。

多くのドメインは、数十〜数百残基程度のアミノ酸からなり(例として、PDZ ドメインは、90残基弱で、ホメオドメインは60残基程度、SAMドメインは70残基程度、SH2ドメインは100残基程度、pairedドメインは100残基以上、MAMドメインは170残基程度...)、固有の立体構造を取り、単体でも何らかの機能を有し、ドメインの組み合わせによって、新規遺伝子が作られることもあったと考えられている。が、これでも、まだちょっと長い。

例えば、１Lの水中で、毎秒1mM相当の新規タンパクが生成する(つまり、毎秒6.02×10^20個の新規タンパクが生成される)ような(法外な)イベントがあったとしても、一億年かけても、最大で10^36個オーダーのタンパク質しか試せない。現在の海水(その量は、10^21リットルのオーダーだそうだ)を全部使っても、20の100乗よりは大分小さく、ランダムな重合によって、100残基のタンパクを全部網羅するには到底及ばない。

というわけで、生物が遺伝子ガチャに成功している理由が説明されるべきで、以下のように、いくつかの可能性が考えられる。
①全ての遺伝子は、単一の遺伝子に由来する
②ランダムに生成したアミノ酸配列でも高確率で何らかの活性を持っている:例えば、ドメインの配列も高度に保存された(連続してるとは限らない)短い配列と、それ以外の部分からなることが多々あり、保存された配列以外の部分は、割と何でもいいとか
③ドメインより小さい構造単位が存在する

現在、多くの遺伝子が、遺伝子重複によって、既存の遺伝子を雛形にして生成し、重複遺伝子は元になった遺伝子と独立な変異を経ることで、しばしば機能を変化させてきたと考えられている。ヒトに多くの相同遺伝子(パラログ)が含まれていることが分かっているし、大腸菌や酵母も、パラログ遺伝子を持っていて、これらの事実は遺伝子重複によって説明できる。

1970年に出版された大野進の本Evolution by gene duplicationには、“In a strict sense, nothing in evolution is created de novo. Each new gene must have arisen from an already existing gene.”と書かれていて、これは①の考え方に相当する。当時は知られてなかった遺伝子の水平伝播も、今ではありふれた現象と考えられてるけど、既存の遺伝子を雛形として新しい遺伝子を獲得するという点では共通している。

進化の歴史上、このような方法によってのみ遺伝子の獲得が起きたとすれば、確率的な困難は、最初の遺伝子が誕生した時に一回だけ存在したことになり、一回だけなら、奇跡が起きることもあるだろう。但し、このような奇跡が必要だったなら、地球以外の星に生命体が存在する可能性は極めて低いことになる。

遺伝子重複は、古くから予想されていたことではあるけど、最近は、個別の遺伝子について、かなり詳細な進化の軌跡を辿ることができるようになっている。

遺伝子重複の顕著な例は、ヒトのLオプシン(赤オプシン)とMオプシン(緑オプシン)で、この2つの遺伝子はX染色体上に並んで存在する。SオプシンはX染色体上になく、分岐したのはもっと古いと思われる(塩基配列上も、最も違いが大きい)。現在、原猿や有胎盤類は、SオプシンとL/Mオプシンの2色性色覚のものが多いらしく、真猿類の中だと、新世界ザルのオスは(一部を除いて)2色型色覚で、旧世界ザル、類人猿、ヒトは3色型色覚を持つ。新世界ザルは、3〜4000年前に真猿類の一部がユーラシア大陸から南米に移住して分岐し、アフリカに残った真猿類の一種で、L/Mオプシンの遺伝子重複が起き、この種は、旧世界ザルと類人猿の共通祖先になったと考えられている。新世界ザルの祖先が、海を渡った方法は謎だけど、当時は今より大陸間は近かったと考えられている。

と、よく書いてあるのだけど、ヒトの場合、X染色体上にオプシンを3つ持つ人がいるらしい。
OPN1LW
https://www.ncbi.nlm.nih.gov/nuccore/NM_020061
OPN1MW
https://www.ncbi.nlm.nih.gov/nuccore/NM_000513
OPN1MW2
https://www.ncbi.nlm.nih.gov/nuccore/NM_001048181
OPN1MW3
https://www.ncbi.nlm.nih.gov/nuccore/NM_001330067

OPN1MW2とOPN1MW3は、目で見ると同じもののようなので、ヒトのX染色体上には、OPN1LWとOPN1MW、OPN1MW2の3つのオプシンがある。OPN1MW2は、ある人とない人がいるらしい(?)から、かなり最近になって出現した遺伝子と考えられる(数千年とか数万年とか)。OPN1MWとOPN1MW2は、アミノ酸配列が一致してるが、一部の人は、これによって4色型色覚となるようなので、多分、なんか多型があるのだろう。

チンパンジーには、(私が探した限り)2つしかないっぽい。ヒトでもチンパンジーでもXq28領域に遺伝子がある。
OPN1LW領域(チンパンジー)
https://www.ncbi.nlm.nih.gov/nuccore/CU467112.1?report=fasta&from=144100&to=156324
OPN1MW領域(チンパンジー)
https://www.ncbi.nlm.nih.gov/nuccore/CU467112.1?report=fasta&from=186659&to=199065
ヒトXq28領域
https://www.ncbi.nlm.nih.gov/nuccore/AC092402.2

どーでもいいけど、
チンパンジーX染色体
https://www.ncbi.nlm.nih.gov/nuccore/1348310337
ゴリラX染色体
https://www.ncbi.nlm.nih.gov/nuccore/NC_044625.1
を見ても、OPN1LW/OPN1MWに相同な遺伝子が一個しか見つからないので、アセンブルをミスってるような気がする。

ヒトのOPN1LWとOPN1MWは、共に364aaで、アミノ酸13個の違いしかない。通常、オプシンは、一つの視細胞に一種類しか発現してないらしいけど、4色型色覚の人で、そのへんの発現制御がどうなってるのかは謎。

OPN1MWとOPN1LWの場合は、上流に共通のエンハンサー領域LCRというのがあって、一方の遺伝子のプロモーターを(ランダムに)選択し、一旦行った選択は変化しないらしい(どうやって安定化するのかは分かってないらしい)。ゼブラフィッシュでは、RH2-LCRという領域が4つの下流にあるRH2遺伝子の発現制御に関与しているそうなので、ヒトでも、OPN1MW2が追加されたところで、同種の機構が働いているのかもしれないけど不明要素が多い

OPN1LW,OPN1MW,OPN1MW2とOPN1SWは、エクソンの数が合致していて、個々のエクソンのサイズもよく似ているから、これらの遺伝子が単一の遺伝子に由来するというのは尤もらしい。ヒトのオプシンファミリー遺伝子には、他にも、OPN2/RHOやメラノプシン(OPN4)などの複数の遺伝子があるが、これらではエクソン数が違っている。

OPN1以外のオプシンファミリー遺伝子も、元を辿れば、単一の遺伝子に起源があると信じられている。更に、オプシンは、GPCR(Gタンパク共役型受容体)で、最初のオプシンは、他のGPCR遺伝子から遺伝子重複と変異による新機能獲得(neofunctionalization)で生成された遺伝子かもしれない。ヒトには1000近いGPCR遺伝子があるが、277個のGPCRの系統樹を作成した論文(PMID 12429062)によれば、オプシンファミリーに最も近いのはメラトニン受容体のようである。カイメンには、メラトニン受容体もオプシン遺伝子もないらしいが、また別の論文(PMID23112152)には、平板動物はオプシンを持つと書いてある。

原核生物や古細菌も、バクテリオロドプシンやSchizorhodopsin、Heliorhodopsinと呼ばれる遺伝子を持つが、これらはGPCRではないようだ。レチナールと結合し発色団としている点は共通し、また7回膜貫通型受容体ではあるらしいから、全てのGPCRの起源という可能性はなくもない(?)(7回膜貫通型受容体という用語は、GPCRとほぼ同義で使われる)

他の例として、ホタルのルシフェラーゼ遺伝子は、遺伝子重複とneofunctionalizationで出現し、ついでに、その後、更に遺伝子重複が起きたと考えられるらしい。
cf)発光生物学　2つの最新トピックス
https://doi.org/10.1271/kagakutoseibutsu.57.409

単一遺伝子の重複だけでなく、ゲノム単位の重複も起きると考えられている。Hox gene clusterは、ショウジョウバエでは、8遺伝子からなる単一のクラスターだけど、多くの脊椎動物では、4つのクラスターがあり、それぞれのクラスターは10個前後の遺伝子からなる(進化の中で、一部が脱落したり増えたりしてるけど、並び順は保存されてる)。ヒトやマウスではHoxa,Hoxb,Hoxc,Hoxdと名前がついてて、異なる染色体上(ヒトではそれぞれ7,17,12,2番染色体)にあるけど、これはゲノム単位の重複が2回起きて生じたという仮説がある。クラスター単位で重複して転座しただけという可能性もあるけど、ParaHoxクラスターとかも4つあるようだ。

他にも、哺乳類のFOXP遺伝子は4つあるけど、ハエには一つしかない。2回ゲノム重複があったと考えると、丁度いいケースが結構ある。2010年代後半までの報告では、脊椎動物に近縁のナメクジウオやホヤが単一のHoxクラスターしか持たない一方、無顎類の生き残りである円口類のHox gene clusterやParaHox gene clusterは複数あるものの、どっちも2つとか、両方4つずつというわけではなく、ゲノム重複の起きた時期は、まだ完全には確定されていない。
cf) Chromosome evolution at the origin of the ancestral vertebrate genome
https://www.biorxiv.org/content/10.1101/253104v1
cf) Lampreys, the jawless vertebrates, contain only two ParaHox gene clusters
https://www.pnas.org/content/114/34/9146

ゲノム重複は魚類や両生類でも起こると考えられ、有性生殖を行うフナ類や金魚は、52対104本の染色体を持ち、1400万年前にゲノム重複があったと2019年に報告された(他の事例は、アフリカツメガエル)。また、多くの真骨魚類が、最低8つのHox gene clusterを持つらしく、数億年前に、初期の真骨魚類でゲノム重複が起きた可能性も考えられている。節足動物では、一部のクモやサソリに2つのHox gene clusterが見られるようである。

哺乳類だと、仮に4倍体が生じても、子孫を残せる可能性は絶望的と思われてるけど、21世紀になって、哺乳類でも、現在知られる中で最大の染色体数を誇るメンドサビスカーチャネズミ(学名:Tympanoctomys barrerae)及び近縁のゴールデンビスカーチャネズミ(学名:Pipanacoctomys aureus)というマイナーな齧歯類が、4倍体の可能性が議論されている。
Discovery of tetraploidy in a mammal
https://www.nature.com/articles/43815

Molecular cytogenetics discards polyploidy in mammals
https://doi.org/10.1016/j.ygeno.2004.12.004

Evolution of the Largest Mammalian Genome
https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5569995/

ヒトの持つ遺伝子の多くが、重複遺伝子として生まれて進化してきたことは確からしいが、全部の遺伝子がそうだと証明されてるわけではない。21世紀になってから、既存遺伝子を雛形としない遺伝子(de novo gene)の生成があるという議論がされるようになった。orphan geneという言葉もあるけど、orphan geneは、遺伝子の水平伝播で獲得した遺伝子でも、そう見えることがあるのに対して、de novo geneは、DNAのnon-coding sequenceから生成されたものと定義される。

de novo gene birthという名前が付いてるけど、ヒト遺伝子の中で、de novo geneだと特定されてるものは今の所ないし、de novo gene birthのメカニズムもよく分かってない。

de novo geneに関する初期の論文として、2006〜7年頃の以下の論文が挙げられる。

Evidence for de novo evolution of testis-expressed genes in the Drosophila yakuba/Drosophila erecta clade
https://doi.org/10.1534/genetics.106.069245

論文中に、We refer to this class of orphan genes as “de novo,” to suggest the possibility that they may derive from ancestrally noncoding sequence.と書かれていて、de novoという接頭辞を使おうと提案したのは彼らっぽい。

彼らは、分子生物学でよく使われるモデル生物Drosophila melanogasterの近縁種Drosophilia yakubaの精巣で発現しているmRNAを解析し、D. melanogaster, D. erecta, D. ananassaeなどの近縁種に相同遺伝子の存在が確認できないD. yakuba固有の遺伝子を探索した。これらの種は、Drosophila melanogaster species subgroupもしくは、Drosophila melanogaster species groupに分類されている。species subgroupとかspecies groupはハエ業界以外ではあんまり使われてないっぽいけれど、属と種の中間の分類階級らしく、とりあえず、Drosophilia melanogasterに最も近縁な種たちを調べたということ。

Wikipediaによると、Drosophilia yakubaに最も近縁なのが、D. santomeaで、次が、D. teissieriとなっている。
Drosophila melanogaster species subgroup
https://en.wikipedia.org/wiki/Drosophila_melanogaster_species_subgroup

こうして見つかった遺伝子は、2〜4個の比較的少数のエクソンからなると予測されている。論文では、これらの遺伝子が、本当に、non-coding sequenceから出現したものであることが証明されてるわけではないし、タンパク質に翻訳されてない可能性もあるが、これらの遺伝子がタンパク質をコードしているならば、そのサイズは、100aa前後で、かなり小さく、突然、巨大な遺伝子が出現する可能性は低いという直感とは整合的である。

2012年の以下の論文では、酵母のORFの内、遺伝子としてアノテーションされてないものを調べている。
Proto-genes and de novo gene birth
https://www.nature.com/articles/nature11184

ゲノム配列からのORFのアノテーションでは、最初に、300塩基を閾値として取捨選択されていて、約6000の遺伝子に対して、約261000のアノテーションされてないORFがあったそうである(過去形だけど、多分、今も状況はそんなに変ってない)(個別に機能が同定された遺伝子の中には、閾値より小さいものもあると思われる)。論文では、長さが30塩基以上のアノテーションされてないORFについてタンパク質に翻訳されてるかどうかを調べた所、108000のORFの内、1139個について翻訳されている証拠を得たとしている("We found that 1,139 of 108,000 ORFs0 show such evidence of translation")。

割合としては1%だけど、そもそもプロモーターがなければ転写もされないので、少ないとも言い切れない(どれくらい転写されてるのかは分からないが)。実験方法としては、リボソームプロファイリングと呼ばれてる方法で、2009年に提案され、次世代シーケンサーを使うので、20世紀には出来なかった実験。但し、翻訳されているというだけで、それが何かの機能を果たしているかは分からない。

"weakly conserved annotated ORF"とアノテーションなしのORFでは、翻訳が確認されるものが1256個、300塩基以下の長さのものは1133個ある(Fig4a及び、aupplementary information fig8 middle)らしいので、100アミノ酸以下の割合が、かなり高い。

2018年の以下の論文は、ちょっと毛色が違って、以前から知られていた遺伝子がde novo geneであると主張している。
De Novo Gene Evolution of Antifreeze Glycoproteins in Codfishes Revealed by Whole Genome Sequence Data
https://dx.doi.org/10.1093%2Fmolbev%2Fmsx311

不凍タンパク質と呼ばれる凍結防止の機能を持つ遺伝子は、様々な生物で見つかっていて、独立に何度も獲得されたと考えられている。論文では、タイセイヨウダラ(Gadus morhua)やコダラ(Melanogrammus aeglefinus)が持つ不凍タンパク質AFGPが、1300〜1800万年前にde novo gene birthで獲得したと推論している。2019年の論文

Molecular mechanism and history of non-sense to sense evolution of antifreeze glycoprotein gene in northern gadids
https://doi.org/10.1073/pnas.1817138116

も、同様の結論だけど、こっちの方が読みやすい(結論だけ見たければ,Fig4)。Fig2(Boreogadus saidaという種のafgp遺伝子)を見ると、3つのエクソンからなり、冒頭2つのエクソンはシグナルペプチドで、3つ目の最も大きな本体は、TAAが繰り返し並ぶ特異な配列をしている。afgp遺伝子は、そこそこ大きいけど、例外的ケースじゃないかと思う。

2019年の以下の論文は、現存する自然界の遺伝子の話ではなく、ランダムに生成した人工的なsmall ORF(長さは10~50aa)によって、大腸菌が抗生物質耐性を獲得しうるという実験。
De Novo Emergence of Peptides That Confer Antibiotic Resistance
https://doi.org/10.1128/mbio.00837-19
10の9乗に満たない種類の配列から、3つの候補遺伝子が得られたそうだ。あくまで、抗生物質耐性の有無でスクリーニングしているだけなので、大腸菌にとって有用な機能を果たす遺伝子は、もっと多く含まれている可能性もある。

ヒトに、de novo geneが、どれくらいあるかは分からないが、現在のところ、小さい遺伝子は見過ごされてる可能性も高いので、意外と多い可能性もある。酵母の例を見ると、既知の遺伝子6000に対して、アノテーションされてないが翻訳されるORFが1000以上あり、これらの内、どれくらいが何らかの機能を果たしているのか分からないが、全遺伝子の一割くらい、de novo geneという可能性もある。

ドメインより小さい構造単位が存在するかどうかについて、以下のような論文を書いてる人たちがいる。

Centripetal modules and ancient introns
https://doi.org/10.1016/S0378-1119(99)00292-9

タンパク質の3次元構造に基づいて、各アミノ酸残基に対して定義される"F-value"(とは呼んでないが、勝手に命名した。周辺アミノ酸のα炭素との距離の二乗の平均)の極小点を境界として、centripetal moduleというものを定義して、moduleの境界は、フェーズ0のエクソン・イントロン境界(エクソンの始まりや終わりが、コドンの途中になってないような境界)と、よく相関すると主張する論文。

centripetal moduleの概念は、1987年に提案されたものらしいが、この論文では、module境界を計算するプログラムを書いて、統計的に相関があるのを確認したということらしい。プログラムの詳細は説明されてないし、公開されてもいない。

彼らは、ドメインについては特に何も言ってないけど、centripetal moduleは、直感的には、空間的にコンパクトにまとまっている単位に分割しようという発想で、ドメインも、しばしばコンパクトにまとまっているとされているので、ドメイン境界を含むものになるだろうと予想される。

色々疑問はあるけど、例えば、ヒトのピルビン酸キナーゼM2(PKM2)で、エクソン・イントロン境界とドメイン境界に、何か関わりがあるか見てみる。

PKM2は、グルコースを2つのピルビン酸に分解する解糖系に於いて、ホスホエノールピルビン酸からADPにリン酸を転移し、ピルビン酸を生成する過程を触媒する酵素である。ヒトは、複数のピルビン酸キナーゼを持つが、その内の一つ。

PKM2(PKM isoform a) domains
http://atlasgeneticsoncology.org/Genes/GC_PKM.html

によれば、PKM2は、4つのドメインからなり、A-domainは、配列上は2つに分断されていて、N-domain (aa 1-43), the A-domain (aa 44-116 and 219-389), the B-domain (aa 117-218) and the C-domain (aa 390-531)らしい。

pkm2遺伝子のエクソンは、
Homo sapiens pyruvate kinase M1/2 (PKM), transcript variant 1, mRNA
https://www.ncbi.nlm.nih.gov/nuccore/NM_002654.6
によれば、11個のエクソンからなるようで、その境界は、タンパク質のドメイン境界と、特に関係してそうには見えない。

PKM2は、解糖系に関わる酵素なので、その進化的起源は、物凄く古い可能性がある。

ヒトPKM2の構造は、
Functional cross-talk between allosteric effects of activating and inhibiting ligands underlies PKM2 regulation
https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6636998/
のFig1(A)で見られ、

Evolutionary plasticity in the allosteric regulator-binding site of pyruvate kinase isoform PykA from Pseudomonas aeruginosa
https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6802521/
のFig3Bに緑膿菌のピルビン酸キナーゼPykAの構造があるが、目で見る限り、ヒトPKM2に似ている。緑膿菌PykAのA-domainには、高度に保存された配列MVARGDLGVEがあり、触媒活性に重要な役割を果たすと書いてある。ヒトのPKM2のA-domainには、似た配列MVARGDLGIEが含まれている(触媒活性に関与するのかは知らない)

それはそれとして、PKM2のcentripetal moduleを調べることにする。構造データとして、
https://www.rcsb.org/structure/1T5A
を使用して、計算した。横軸は、残基で、縦軸は、F-value。kは、周辺アミノ酸何残基までを計算に入れるか決めるウインドウ幅(の半分。論文では、k=40~90がいいとか書いてるので、同じような範囲で計算した)。
f:id:m-a-o:20201213154649p:plain
は、平滑化なしの値をプロットしたもので、既知のドメイン境界も、重ねて書いてある。

PKM2のドメインは、目で見ても、何となく分かるので当然のことではあるが、ドメイン境界が極小点になってそうだとは分かる。この例では、極小点の中でも、ドメイン境界になっているものは、F-valueがウインドウ幅kの選び方に依存してないようである。いつでも、こんな風にドメイン境界が採れるという保証はない。

生データだと、極小点が多すぎるので、8残基移動平均をとって平滑化して、極小点を抽出してみた。この例だと、あまり、エクソン・イントロン境界と相関してそうに見えない。"+"が付いてるのは、フェーズ0で、"-"が付いてるのは、フェーズ0じゃないもの。
f:id:m-a-o:20201213154709p:plain
平滑化したグラフに、ドメイン境界を重ねてみる
f:id:m-a-o:20201213154712p:plain

もうちょっとencourageされる結果だったら、他の例も見てみようと思ったけど、大変微妙な結果。今の所、10〜20残基程度の短い構成単位が存在するという証拠はなさそうに思える。

描画に使用したコード

import matplotlib.pyplot as plt
import numpy as np
from scipy import signal

def readcif(filepath):
   ret = []
   nprev = None
   for line in open(filepath):
       line = line.strip()
       if len(line)<5:continue
       if line[:5]!="ATOM ":continue
       ls = line.strip().split()
       if len(ls)<10:continue
       if ls[0]!="ATOM":continue
       if ls[3]!="CA":continue  #-- not alpha carbon
       if ls[4]!="." and ls[4]!="A":continue
       if ls[6]!="." and ls[6]!="A":continue
       aa = ls[5]
       n = int(ls[8])
       x,y,z=float(ls[10]),float(ls[11]),float(ls[12])
       assert(nprev==None or nprev+1==n),(nprev,n)
       ret.append( (aa,x,y,z) )
       nprev=n
   return ret


def getseq(r):
   dic = {'Ala': 'A', 'Asx': 'B', 'Cys': 'C', 'Asp': 'D', 'Glu': 'E', 'Phe': 'F', 'Gly': 'G', 'His': 'H', 'Ile': 'I', 'Lys': 'K', 'Leu': 'L', 'Met': 'M', 'Asn': 'N', 'Pro': 'P', 'Gln': 'Q', 'Arg': 'R', 'Ser': 'S', 'Thr': 'T', 'Sec': 'U', 'Val': 'V', 'Trp': 'W', 'Xaa': 'X', 'Tyr': 'Y', 'Glx': 'Z'}
   def normalize(s):
      s= s.lower()
      return s[0].upper()+s[1:]
   return "".join([dic[normalize(x[0])] for x in r])



def getF(r , k=50):
   ret = []
   for i,(aa,x,y,z) in enumerate(r):
       j1 = max(0 , i-k)
       j2 = min(len(r)-1 , i+k) + 1
       assert(j2>j1),(j1,j2)
       Fval = sum([(x-x2)**2+(y-y2)**2+(z-z2)**2 for(_,x2,y2,z2) in r[j1:j2]])/(j2-j1)
       ret.append(Fval)
   return ret


def phase(n):
  if n%3==0:return "+"
  else:return "-"


if __name__=="__main__":
  r = readcif("1t5a.cif")
  Nlen = 531  #-- PKM2 has 531aa
  offset = Nlen - len(r) + 1
  xtick = "module"
  smoothing = True
  for k,c in [(40,'r'),(50,'brown'),(60,'grey'),(70,'g'),(80,'b')]:
    d = getF(r,k=k)
    x = np.linspace(offset ,Nlen, len(d))
    y = np.array(d)
    if not smoothing:
      plt.plot(x,y,c,label='k={0}'.format(k))
    else:
      num=8  #移動平均の個数
      b=np.ones(num)/num
      y2=np.convolve(y, b, mode='same')#移動平均
      plt.plot(x,y2,c,label='k={0}'.format(k))
      minpos = signal.argrelmin(y2, order=3)
      plt.plot(x[minpos],y2[minpos],c, marker='o',linestyle='None')
    if xtick=="exon":
       exon_ticks=[((c-89)/3,str(n+3)+phase(c-89)) for (n,c) in enumerate([243,335,467,654,925,1076,1229,1396,1578])]
       plt.xticks([c[0] for c in exon_ticks],[c[1] for c in exon_ticks])
    elif xtick=="module":
       plt.xticks([44,117,219,390])
    plt.yticks([])
    plt.grid(True)
    plt.xlim(0,Nlen)
  plt.legend()
  plt.show()