1Dp01
電子移動を利用した DNA 塩基配列検出のモデル(東芝・豊橋技科大・科技団・筑波大)○田中成典・夏目貴行・薄木亮・仙石康雄・
小川哲司・栗田典之・横島智・岡田朗
DNA 鎖を分子ワイヤーとして用いる応用の一つとして、DNA における一塩基多型（Single Nucleotide Polymorphisms；SNPs）を検出するための電気化学的 DNA チップがある。このデバイスにおいては、従来の蛍光方式の DNA チップと異なり、一本鎖と二本鎖の違い、あるいは SNPs による塩基ミスマッチの存在を電気化学的応答の差として検出する。本研究では、DNA の構造や塩基配列の違いによる電子移動速度定数の差を記述する理論モデルを構築し、電気化学的 DNA チップに対する実験結果との定量的な比較を試みる。
カリフォルニア工科大学の Barton らのグループは、図１に示すような金電極上に DNA
鎖を固定した構造を用いて塩基ミスマッチの電気化学的応答に対する影響を系統的に測定した [1]。彼らの結果によると、AC ミスマッチの存在によって相補的な AT ペアの場合と比較して移動電荷量の大幅な減少が見られる（図２）。この実験結果を説明するために、我々は図３に示すようなモデルを考え、DNA 鎖を介した電極とインタカレータ間の電荷移動速
度定数と移動電荷量を計算した。電極電位に依存して電子移動メカニズムは直接的（コヒーレント）な超交換機構から間接的（インコヒーレント）なホッピング機構へのクロスオーバー（図４）を示し、図５に示すように、AC ミスマッチの存在によって電荷移動速度定数の低下が見られる [2]。但し、この解析においては、AC ミスマッチの存在によってその周囲の電子移動積分が約半分程度に減少するという現象論的な仮定がなされた。
本研究では、より第一原理的な立場から実験との比較を行うために、まず分子力学および分子動力学法によって水溶液中の DNA 鎖の構造を求め、それを用いて分子軌道法により電子状態計算を行い、電子移動反応速度を評価するために必要な分子パラメータを算出するという方法を取った。DNA 鎖の構造計算は相補的な塩基対ならびにミスマッチのある場
合（図６）の双方に対して行われ、さらに二本鎖と一本鎖の比較も行った（図７）。また、Protein Data Bank 等に登録されている実験構造との比較も行った。
図８は AC ミスマッチ部分を切り出して計算した HOMO（A に局在）、LUMO（C に局在）の電子分布である。分子軌道計算に基づきミスマッチ近傍の電子移動積分を評価してみると、確かに AC ミスマッチの場合に AT 相補対の場合と比べて移動積分の値が低下する
ことが確認され、以前のモデルの正当化がなされた。このように、構造→電子状態→電子移動反応速度の計算の道筋をボトムアップ的に辿ることにより、塩基ミスマッチの存在を
電気化学的に検出する手法の基礎付けを理論的に行うことが可能となりつつある。
本研究は、科学技術振興事業団・計算科学技術活用型特定研究開発推進事業の研究課題「DNA のナノ領域ダイナミクスの第一原理的解析」の援助を受けて行われた。

e

図 1．電気化学的 DNA チップの模式図

図３．DNA チップのモデル

図２．
AC ミスマッチによる電荷移動量の低下 [1]図４．
電子移動スキーム

図５．速度定数の計算結果 [2]Radius of Inertia (A)

13.5
12BP Double Strand
13.0
12.5
12.0
12BP Single Strand
11.5

図６．ミスマッチがある場合の MD

0

100
200
Time (ps)300

図７．２本鎖と１本鎖の慣性半径の時間変化

図８．AC ミスマッチ塩基対の HOMO（左）ならびに LUMO（右）［文献］[1] S.O. Kelley et al., Nucleic Acids Res. 27 (1999) 4830.
[2] A. Okada et al., J. Mol. Struct. (Theochem) 630 (2003) 283.

1Dp02ヘムオキシゲナーゼによるヘムの代‫ޥ‬機構に関する理論的研究(九大先導研、Institute of Problems of Chemical Physics†)              ○蒲池݄志、Alexander F. Shestakov†、吉澤一成【序】ヘムオキシゲナーゼはヘムをビリベルジン、CO、鉄とに分ӂする酸素添加ܼ素である。動物では鉄のリサイクル、植物では光センサーに重要なビリン໸の前駆体であるビリベルジンの生成などの役割を果たしている。Scheme に示すように、本ܼ素によるヘム分ӂ過程の中間体にはa-ヒドロキシヘム、ベルドヘムが知られており、a位の炭素が選択的に酸化される。このようにヘムオキシゲナーゼは酸素添加反応を触媒するにも関わらず単純タンパク࠽であり、ヘムは基࠽であると同時に酸素活性化のための補欠分子族として働く。P450 では反応活性種は鉄オキソ種(Fe(IV)=O)
であるが、ヘムオキシゲナーゼの反応ではヒドロパーオキソ種(Fe(III)–OOH)が反応に関わると考えられてきた。‫ة‬年、Hoffman 等は EPR、ENDOR の測定からこの反応機構を支持している[1]。本研究では、この反応機構を密度汎関数法を用いてӂ析を行った。
Scheme
Me

V

a

V
N
Fe

d

Me

Me
V
NADPH/O2

N

N
g

Me

V
Me

N

N

V
O2

H2O

Pr

heme

N

N

Me

N

V
OO

Me
V
NADPH/O2

Me

NH HN

Fe
CO

N

Me

Me
Pr

V

+
O

Fe

b

N
Me

Pr

OH

Me

Pr

a- meso -hydroxyheme

Me

N

N

N
Me

Pr

Pr

verdoheme

HN

Me

Me
Pr

Pr

biliverdin
V = vinyl Pr = propionate【‫ٽ‬算方法】‫ٽ‬算方法として密度汎関数法のひとつである B3LYP 法を用いた。基底関数には全ての原子に 6-31G を用い、さらに Fe、O、N、a-carbon には d 型分極を加えた。࠰配位子のモデルとしてイミダゾールを採用した。ヒドロパーオキソ種の‫ٽ‬
算にはஏ荷を 0、スピン多重度を基底状態である 2 重項とし、オキソ種の‫ٽ‬算にはஏ
荷を 1、スピン他重度には２重項と４重項のふたつを考慮した。全ての‫ٽ‬算にはGaussian98 を用いた。
【結果と考察】図 1 に‫ٽ‬算で得られたヘムオキシゲナーゼのヒドロパーオキソ種によるヘムのa-ヒドロキシヘムへの転換反応のエネルギーダイアグラムを示す。ヒド
ロパーオキソ種の最適化構造(R-hp)の O–O 結合ଥは 1.446 Å であり、CH3S-が配位しているシトクロム P450 のヒドロパーオキソ種のモデル錯体に比べ 0.094 Å 短いものであった。この O–O 結合ଥの変化は Dawson 等が提案したシトクロム P450 の࠰配位子であるシステインの強いஏ子供与

2.147

1.940

1.630

性のためと考えられる。

2.217

反応の初期段階において、ヒドロr (Fe–N)= 1.982

パーオキソശの OH 基がヘムのa-炭

1.792
1.050

TS1-hp42.9

1.728

素に直接攻撃し、中間体 I1-hp が生

1.392

1.525
1.381
2.072

1.437
1.186
r (Fe–N) = 2.021

2.103
r (Fe–N) =1.991

成する。その後、転移した OH 基のTS2-hp8.3水素が生成したオキソによって引き0.0

抜かれる遷移状態 TS2-hp を経て、

7.3
1.796

1.446
1.785

–10.7

2.103

1.430
1.627

ィリン上の酸素がa-炭素から水素を

1.334

2.067

1.920

I2-hp が生成する。さらに、ポルフ

TS3-hp18.8

r (Fe–N) = 2.013

I2-hp

2.179
r (Fe–N) = 2.018

引き抜くことでa-ヒドロキシヘムが

R-hp

–48.4
r (Fe–N) = 1.998

P-hp

I1-hp

生成する。
また、シトクロム P450 の反応で主な活性種として働く、オキソ種についてもӂ析を行った。図 2 に示す
ように、ポルフィリンが大きく歪む

2.080
r (Fe–N) = 2.020

図 1. ヘムオキシゲナーゼのヒドロパーオキソ種によるヘム酸化過程に関するエネルギーダイアグラ
ム。(単位は kcal/mol、Å)2.042 (1.946)

ことで、オキソがa-炭素と結合

1.697 (1.711)

し、中間体 I1-oxo が生成する。

2.116 (2.111)

この歪みを緩和するため、Fe–O

r (Fe–N) = 1.971 (1.966)

42.5

1.404 (1.437)
1.208 (1.194)

r (Fe–N) = 2.005 (1.989)

TS2-oxoTS3-oxo

34.5

39.9

Doublet35.3

Quartet
19.3

23.7

16.4
1.438 (1.408)‫ٽ‬算の結果得られた活性化エ
0.0

オキソ種ともに 40 kcal/mol ほど

1.377 (1.377)
2.030 (2.129)

TS1-oxo

生成物に至る。
ネルギーはヒドロパーオキソ種、1.405 (1.399)

2.805 (2.507)

r (Fe–N) = 1.967 (1.961) 2.166 (2.256)

結合が開裂する。その後、ヒドロパーオキソ種の場合と同様に

1.370

1.793

19.6

25.9

14.8

1.890 (2.041)
1.367 (1.351)–0.2

2.043 (2.409)

2.170 (2.132)

1.619 (1.619)
r (Fe–N) = 2.007 (1.981)
2.163 (2.161)と極めて݄い。例えばシトクロr (Fe–N) = 1.955 (1.976)

I2-oxoI1-oxo–36.3

r (Fe–N) = 2.024 (2.024)

ム P450 の Compound I による強

R-oxo
–45.8

固な C–H 結合の引き抜き反応で1.928 (2.131)

さえ半分の 20 kcal/mol ほどである[2]。本研究の結果は従来考えられてきた反応機構について再1.353 (1.356)

r (Fe–N) = 1.991 (1.997)

P-oxo

図 2. ヘムオキシゲナーゼのオキソ種によるヘム酸化
過 程 に関 す る エネ ル ギ ーダ イ ア グラ ム 。 ( 単 位 はkcal/mol、Å)考を促すものとなった。(J. Am.
Chem. Soc. submitted.)[1] Davydov, R.; Kofman, V.; Fujii, H.; Yoshida, T.; Ikeda-Saito, M.; Hoffman, B. M. J. Am. Chem. Soc. 2002, 124, 1798. [2]Kamachi, T.; Yoshizawa, K. J. Am. Chem. Soc. 2003, 125, 4652.

1Dp03Ebtellur 誘導体の Thiol Peroxidase 反応の機構解明（東大院工 A・Emory 大 B）

○崎元

友紀子 A、Djamaladdin G. MusaevB、平尾 公彦 A

【序論】Peroxynitrite（PN）は生体内における強力な酸化剤であり、DNA 損傷や脂質の過酸化の誘発、酵素の不活性化など、人体に重大な影響を与える。これに対し、様々な有機 Te化合物が有機 Se 化合物と同様に抗酸化剤として働き、Thiol

peroxidase 活性を示すことが知られている。触媒反応のスキームを図１に示す。触媒サイクルは PN による Te の酸化反応（第一段階）とチオールによる Te の還元反応（第二段階）との大きく 2 つに分かれる。
-(第一段階)

EbTe + ONOO -EbTe=O + NO 2

(第二段階)

EbTe=O + 2RSH

EbTe + RSSR + H 2 O

しかしながら、その反応メカニズムの詳細については未だわかっておらず、これらの解明が高活性な触媒設計に不可欠である。Ebtellur 誘導体(図 2)はチオールを用いて PN を還元する
触媒のひとつであり、その擬物質 Ebselen は抗炎症剤として注目を集めている。また Te 化
合物はその Se 化合物よりも高活性を示すものが多いという報告があり、大変興味深い。そこ
で、本研究は Ebtellur 誘導体と PN の触媒反応についての反応機構を明らかにすることを目的としている。
RSSR
+
H2O

EbTe

O
ONOO -

N R

RSH
ONO

-

M

OH
EbTe

EbTe=OSR
RSH

図 1 反応スキーム

図2

Ebtellur 誘導体【計算方法】分子軌道計算には密度汎関数法（B3LYP）を用い、基底関数系は Te に対しては Stuttgart-Dresden relativistic ECP を、その他の原子に対しては 6-311+G(d,p)を用いて、構造最適化及び振動解析を行った。チオールには CH3SH を用いて計算を行った。
【結果】（第一段階）図３に第一段階の反応のポテンシャルエネルギーを示す。反応は EbTe +PN → EbTe-PN → EbTe-TS1 (O-O activation) → EbTe-O(NO2–)のように進行することがわかった。Ebtellur 誘導体と PN が反応し、EbTe-PN の錯体が形成された後、O-O 結合が活性化され、遷移状態 EbTe-TS1 を経て、中間体 EbTe-O(NO2–)が形成される。活性障壁は 9.8kcal/mol
であった。中間体 EbTe-O(NO2–)から反応経路は NO2–の解離と NO3–への異性化との二つに別れた。NO3–への異性化は遷移状態 EbTe-TS２(NO3- formation)を経て進行し、中間体 EbTe-NO3–が生じた後、NO3–の解離が起こる。NO3–の異性化の活性障壁は 26.9kcal/mol であった。NO2–の解離エネルギーは 25.9kcal/mol であり、両方の経路とも EbTe-O(NO2–)からの段階、それぞれ NO2–の解離と NO3–形成が律速段階であった。NO2–の解離エネルギーと NO3–の異性化の活性障壁との差は 1.0kcal/mol とわずかなので、2 つの反応が同時に起こると考えられる。また、反応において Te-N 結合は非常に柔軟性に富んでいた。弱い Te-N 結合は PN の錯形成、O-O
結合決裂、NO3–-異性化に有利であり、強い Te-N 結合は NO2–への解離に有利であると考えられる[1]。
∆H (∆G)
0.0

-10.0

EbTe
O

-20.0

Structures

EbTe-O(NO 2-)

EbTe-TS1

EbTe-PN

EbTe-TS2
EbTe -O + NO 2-

O

ONO

N

EbTe

EbTe

O

O

O

-29.2 (-18.8)

NH

OONO

-30.0

-36.2 (-25.0)-39.9 (-28.6)+ NO 2 -Te-26.2 (-26.0)

O

-40.0

O
NH
-62.8 (-51.9)

-50.0

Te

ONO

-52.7 (-52.2)

O

EbTe + NO 3 -

EbTe

-60.0

+ NO 3 --79.3 (-68.1)

-70.0

O N

EbTe

O
O

EbTe-NO 3-

図 3 第一段階の反応のポテンシャルエネルギー(第二段階）図 4 に第二段階の反応のダイアグラムを示す。第一段階の NO2-の解離後生じた、
Telluroxide(EbTe＝O)にチオールが挿入され、遷移状態を経て中間体１Te(OH)(SCH3)が生成する。次にこの中間体１と第二のチオールとが反応し、Disulfide の還元的除去により、EbTe=O
は Ebtellur 誘導体へと還元される。第二のチオールとの反応についてはチオールの S が直接中間体１の S を求核攻撃する経路と Te を求核攻撃する経路との２つの反応経路に分岐する。
計算結果の詳細については当日報告する。
O
N H + RSHTe
EbTe=O

EbTe

EbTe SR

O HSR

O

H

EbTe SR

EbTe

SR

O

SR

O

HO H SR

H

H

OH
EbTe
SR
中間体1
EbTe

SR

SR

EbTe
SR

EbTe

SR

SR
- RS-SROH+ RSH EbTe SR
EbTe
HO
SR
HSR
中間体1

SR

EbTe SRHO H SRNH
Te

EbTe SRHO H

O- RS-SRSR

図 4 第二段階の反応のダイアグラム

R=CH3[1] Y. Sakimoto, K.Hirao, and D.G.Musaev, J. Phys. Chem. A 107, 5631 (2003)

Ebtellurderivative

1Dp04

SAC-CI 法による生体分子のスペクトロスコピー:
Green Fluorescent Protein と Phytochrome(京大院工) ○長谷川

淳也・藤本

和宏・中辻

博

1. Introduction
光が関わる生体内の主要な化学反応過程において、分子構造の変化に伴い光吸収スペクトルが変化する例は数多く知られている。このような系では、構造と励起スペクトルの関係を直接的に研究できる理論計算の重要性は高い。我々はこれまで SAC-CI 法を用いて生体分子の光吸収スペ
クトルの研究を行い、同方法から信頼性の高い帰属が得られることを示してきた。最近の SAC-CI
法のプログラム開発の結果、より大規模な系での計算時間が短縮し、応用性を高めることができた。実際に、対称性の無い生体分子系において active orbital が約 450 軌道程度の計算を実行している。本発表では、Green Fluorescent Protein (GFP)と Phytochrome における光化学に SAC-CI 法を応用し、それらの光化学過程における分子構造と光吸収スペクトルについての研究を報告する。
2. Structure and excited states of GFP chromophore and its mutants
GFP はオワンクラゲ(Aequorea victorea)に含まれる緑色蛍光タンパク質であり、その優れた蛍光特性により分子生物学における分子マーカーとして広く用いられている。また、吸収・発光波長の異なる幾つかの mutant が作られている。これらの光化学を研究する上では、基底状態の分子構
造やタンパク質との電子移動過程を理解することが必要となっている。
Fig. 2 と Table 1 に GFP と mutant の構造と励起・蛍光エネルギーを示す。SAC-CI 法の結果は実験値をよく再現している。これまで BFP と CFP の基底状態の構造が不確定であったが、今回の結果は、それぞれ BFP-II と trans-CFP の構造をとることが示唆している。
また、Thor らにより GFP は UV 光(4.9eV)の照射により、GFPUV へと非可逆的に構造変化することが示され、Glu222 から色素への電子移動が引き金にとなる仮説が提案されている。そこで光化学過程に関与するとされるアミノ酸残基を全て含めた計算モデルを構築し(Fig. 3)、SAC-CI 法により励起状態を求めた。基底関数は DZP に anion 基底を含め、active space は 447 次元とした。Table
2 に示したように、励起エネルギーで 5.6eV までの 10 状態が算出され、51A1 状態（4.18eV)が Glu222にも振幅を持つ軌道からの励起である。しかし、population 解析の結果、移動した電荷はごく僅かであり、GFP 色素から Glu222 への CT 状態の性質は小さいことが分かった。その他の励起状態は色素に局所的な励起状態であった。即ち Thor が提唱する仮説が否定される結果が得られた。
(a)GFP(b)BFP-IO

Table 1. GFP とその変異体の色素の励起・蛍光エネルギーO

N CH3 HNHON(c)BFP-IINN(e)cis-CFP

CH3
O

HN

N

(f)YFP

CH3

Exptl.

SAC-CI

Exptl.

Proteinb

Protein

Gas

Protein

3.12

2.70

2.70 (緑)SAC-CI

GFP

3.23

3.21

BFP-I

3.41

3.16

CH3

BFP-II

3.11

3.22

cis-CFP

3.15

N CH3
N

CH3

Fig. 2 Structure of GFP chromophore and its mutants.

aFluorescence energy (eV)

Gasa

N CH3

O

N CH3
N

CH3

O

N CH3
NH

Molecule(d)trans-CFP

O

N

HN

N

CH3

Excitation energy (eV)

N CH3

2.92
3.25

2.77

2.78 (青)2.68trans-CF

2.99

2.84

2.55

2.56 (紺)

YFP

3.33

3.44

2.81

2.81 (青)Gas phase calculations. bElectrostatic effect is included by point-chage model.

Table 2

Excited states of the active center of GFP.
SAC-CI

State
Character

Exptl.
Eex (eV)aOsc. (au)21A1

GFP π→GFP π*

3.27

0.56

31A1

GFP π→GFP Ryd.

3.98

0.40×10-2

1

4 A1

GFP π→GFP Ryd.

4.11

1.7×10-3

51A1

GFP π, Glu222→GFP π*

4.18

2.7×10-2

1

6 A1

GFP π→GFP Ryd.

4.34

1.1×10-2

71A1

GFP π→GFP π*

4.56

1.2×10-2

1

8 A1

GFP π→GFP π*

4.85

0.15

91A1

GFP π→GFP Ryd.

4.95

6.8×10-3

101A

GFP π→GFP Ryd.

5.17

1.0×10-2

111A

GFP π→ GFP Ryd.

5.58

8.9×10-2

aExcitation energy. b Oscillator strength..

cH. Lossau et al. Chem. Phys. 213 (1996) 1-16.

b

GFP Chromophore

Eex (eV)
3.12c

W22
Glu222
Ser205
Fig. 3 Computational model of GFP active site

4. Excited states of linear tetra-pyrroles : Photo-isomerization mechanism of phytochrome

ビリン化合物はポルフィリン類と同様にテトラピロール構造

H H H Hを有し、アンテナ蛋白質の光吸収や光受容蛋白質の光センサーとして、重要な役割を担っている。本研究ではビリン化合物としてアンテナ蛋白質に含まれるフィコシアノビリンと、光受

O

A

H

N
H
H

容蛋白質、フィトクロムに含まれる色素分子フィトクロモビリン
について基底状態、励起状態の電子構造、及び光吸収スペ

H

H

N
H

N

B
H

D

H

H
N

C
H

H

H

Fig. 4 フィコシアノビリンの基本骨格

クトルに関する研究を行った。
まず、フィコシアノビリン(Fig. 4)とクロリン（環状テトラピロール）について分子軌道と励起状態を比較し、開環効果が励起状態に及ぼす影響を検討した。クロリンのスペクトルは Gouterman の 4 軌道モデルにより、Q-、B-帯と呼ばれる特徴的なピークが現れる。一方、フィコシアノビリンの分子軌道は開環効果や末端の酸素原子置換により、4 軌道モデルは完全に失われる。
フィトクロムは赤色・遠赤色光受容蛋白質であり、形態形成など植物の様々な生理反応に深く関わっている。フィトクロムには２つの異性体（Pr, Pfr）が存在する。生理的に活性な Pr は赤色光（666 nm）により、不活性な Pfr に異性化し、逆に Pfr は遠赤色光（730 nm）により Pr に異性化する。これらの変換はフィトク
ロム中の色素分子であるフィトクロモビリンの光異性化が引き金になっている。これまで２種類の光異性化
機構、１）Pr:ZZZasa→Pfr:ZZEass、２）Pr:ZEZaas→Pfr:ZEEaaa が提案されたが、結晶構造は未だ得られておらず結論に至っていない。本研究では各異性体について SAC-CI 法により励起スペクトルを求め、実験結果と比較し、光異性化機構と色素のプロトン化状態ついて検討した。
色素をプロトン化すると励起エネルギーは大幅に
赤方シフトするので、蛋白質中では中性に近い構造をとると考えられる。また、機構２に対応する異性体
は実験の吸収位置から 0.3eV 程度ずれており、Pr とPfr の吸収強度の関係も再現できなかった。機構１に対応する異性体の励起スペクトルが実験結果に近い値を示した(Fig. 5)。従って色素の異性化は機構１によることを示唆する結果が得られた。
Fig. 5．Pr, Pfr の光吸収スペクトル。

O

1Dp05
ＧＲＩＤコンピューティングによるタンパク質の量子化学計算― ＱＭ／ＭＭ連成プログラムの開発と将来展望 ―（阪大蛋白研、阪大院理１、日立東日本２、日立基礎研３、ＮＥＣ基礎研４）
米澤康滋、中村春木、山中秀介１、山口兆１、黒澤隆２、何希倫３、
佐久間俊広４、中田一人４、〇高田俊和４【序】次世代の計算環境として、GRID が注目されている。GRID では、電力の供給と同じように、パソコンをインターネットに接続するだけで、必要なコンピュータパ
ワーの供給が受けられるようになるといわれている。スーパーコンピュータが手元になくても、例えば、タンパク質の計算が、ブラウザとの会話だけで行えることになる。
これは、シミュレーションスタイルの変革であり、分子シミュレーションのみならず、色々な分野のシミュレーションが、遠からずグリッドコンピューティングを指向するものと考えられる。
このような背景の下、大阪大学バイオグリッドプロジェクトでは、GRID 技術のバ
イオインフォマティクスへの応用を目指して、１）グリッド基盤技術、２）データグ
リッド技術、３）コンピューティンググリッド技術、４）遠隔操作技術、の４分野について、緊密な連携を取りながら、研究開発を進めている[1]。コンピューティンググ
リッドグループの研究活動のひとつとして、分子動力学、分子軌道法、密度汎関数法を連成した QM/MM（Quantum Mechanics/Molecular Mechanics）プログラムの開発を進めている。本講演では、そのプログラムの開発方針と構造的特徴、GRID 環境
下での稼動状況、生体分子への応用などについて、報告する。尚、生体分子のシミュ
レーションに必要な連成プログラムを作成するための基本な考え方を、BioPfuga(Biosimulation Platform United on Grid Architecture) と呼んでいる。
【プログラムの構成】ＭＭ計算を行う prestoX-basic と QM 計算を行う AMOSS を連成したプログラムの構造と、制御及びデータの流れを、図に示す。本プログラムでは、
１）AMOSS が QM 空間に属する原子核に働く力を QM・MM 空間の相互作用も考慮して計算する、２）prestoX-basic が MM 領域の原子核に働く力の計算と AMOSS の結果を合わせ、δｔ秒後の新しい原子核の座標を計算する、３）Couple Controller
なるプログラムが、これらのデータの流れとプロセスの生成を、同期を取りながら制御する、のステップを繰り返して、MD 計算を行っている。GDFT プログラムによる電子相関の計算は、QM 計算部に組み込まれて行われるが、現在この連成作業を進めている。この構造の為、prestoX-basic と AMOSS を、GRID 上で同時並行的に稼動させることが可能であり、プロセッサが十分あれば、QM・MM 空間に対する並列処理による高速化が期待される。また、本システムでは、連成を迅速に行うために、データの入出力は、それぞれのプログラムで独立に行っているが、連成が進めば、自然に融合されていくと考えている。
【プログラムの特徴】 GRID では、分散処理が基本である。また、連成プログラムの構成単位間のインターフェースを簡素化するために、
１）中間データに統一された名前をつけ、ファイルに保存することで、中間データの共有化を図る。現時点では、バイナリデータ中心であるが、一部に関しては XMLを用いている。データのリーダビリティの向上には、全てのデータを XML 化す
る必要があると考えている。
２）プログラムを、個別計算機能に分割し、所謂 SPMD から MPMD への移行を行っている。spawn との併用により、ダイナミックなプロセス生成が可能となり、GRID に適したプログラム構造となっている。
の２点を、プログラム開発の基本方針としている。長期間に渡り、独立に開発されてきた prestoX-basic と AMOSS が、極めて短期間に融合されたことにより、引数によるインターフェースに比べて、格段に連成し易くなっていることが示された。
【計算結果と今後の課題】エタノール＋水系の QM/MM 分子動力学計算を、動作確認のためのテストとして行った。QM 空間として、エタノールと２分子の水、ＭＭ空
間には２２６個の水を配置して、分子グラフィックスによりその動きを表示した。現在、GDFT を連成しているので、電子相関を考慮した MD 計算を計画している所である。また、aquaporin など膜タンパク質への応用も開始している。BioPfuga に基づく連成プログラムは、極めて柔軟性に富んでおり、色々なプログラムと連携しながら、発展していく可能性が高いと考えられる。今後，多くの賛同者を得て、生体分子の汎用シミュレーションシステムとして、機能強化が図られることを期待している。

［１］http://www.biogrid.jp

1Dp06
局在化軌道 QCLO によるバラエティに富んだタンパク質全電子計算（東大生研計算科学技術連携研／アドバンスソフト）○柏木浩、（東大生研計算科学技術連携研）上野哲哉、佐藤文俊、（富士総研）村松伸哉、（アドバンスソフト）稲葉享、（富士総研）西川宜孝【序】
私たちのグループは Gauss 型関数をベースとする密度汎関数法プログラム ProteinDF を開発し、これを用いて、シトクロムｃなどタンパク質の全電子計算を行ってきた[1]。この研究の過程
で私たちが提案した擬カノニカル局在化軌道（QCLO）は、Kohn-Sham 方程式の収束を改善するばかりでなく、複数の金属を持つタンパク質の複雑なスピン状態、巨大タンパク質の分割計算、
溶媒の水の近似計算など、バラエティに富んだ応用が可能である[2]。これらについてのテスト計
算の結果を報告するとともに将来の可能性を検討する。
【QCLO[2]】QCLO は局在化軌道（LO）の一種であり、例えばタンパク質を構成するアミノ酸や、その一部分のペプチド鎖、活性中心といった任意の範囲で局在化しているが、その範囲内ではカノニカル
分子軌道の性質を持つ軌道である。QCLO は次のように求める。
計算対象のタンパク質を小さな単位に分割し、それぞれの分子において通常の分子軌道計算を実行する（本研究では ProteinDF を用いている）。次に、それら分子のカノニカル分子軌道を結合や孤立電子対に局在化するような LO に変換する。続いて、側鎖やペプチド結合の部分といった適当な範囲（フラグメントと呼ぶ）の LO を選び出し、Kohn-Sham 行列や Fock 行列を変換して固有値方程式を解くことによって得られた軌道が QCLO である。そのため、QCLO はフラグメ
ントに局在化しているが、カノニカル分子軌道の性質を持ちフラグメント全体に広がっている。
一度、QCLO を求めてしまえば、QCLO を変換行列として使い、固有値方程式をフラグメントごとに解く操作を分子軌道法の SCF 計算内部で行うことで、フラグメントを簡単に拡大することが出来る。
【シトクロム c 活性中心近傍モデルへの適用】
シトクロム c の鉄ポルフィンと上下に第 5・6 配位子を含む各 3 残基ペプチド鎖 Cys-His-Thr,
Lys-Met-Ile が配位した活性中心近傍モデル（図 1）の計算を行った。3 つのフラグメントに分割した QCLO とカノニカル分子軌道との比較をおこなったところ、QCLO はフラグメント内部でカノニカル分子軌道と極めてよく一致することが確かめられた。
この結果から、QCLO が複数の活性中心を持つ系の計算に応用できることが示唆された。すな
わち、複数の活性中心を別々のフラグメントに属させることによって、それぞれの活性中心においてスピン状態を適切に指定できるため、系全体として複雑なスピン状態が計算できる可能性がある。

図 1：シトクロム c 活性中心近傍モデルの構造【インスリンの自動計算】
QCLO をつなげて、タンパク質の高品位な初期値を作成することができる。タンパク質の密度汎関数法による計算をより安全に達成する手段として有効で、同時に SCF 計算の繰り返し回数が低減する。この方法の有効性を示すため、インスリンの全電子計算に適用し、無理なく全電子計
算を達成している。インスリンはジスルフィド結合が複雑に絡まりあったヘテロダイマータンパ
ク質である。valence double 型の基底関数セットを使用し、軌道数、電子密度の補助基底関数の数はそれぞれ 4443、8074 であった。ポイントは、複雑なトポロジーを生むジスルフィド結合に対応するために、ジスルフィド結合を構成する Cys を持つフラグメントでは結合相方の Cys もつ
なげた状態でフラグメントを定義することである。
現在、極端に計算量を増やすことなく、タンパク質の溶媒分子を考慮にいれた計算を行う方法
として、水分子を QCLO として取り扱う方法をテスト中である。
【QCLO 計算の並列処理】
フラグメント内で局在化しているとはいえ、タンパク質のような大型の系では依然多くの計算時間を要する。そのため、QCLO 計算の並列処理による高速化は欠かせない。問題は、上記の例
で見たように、フラグメントの数やサイズは、活性中心近傍から水分子までまちまちであり、単純な並列処理では高効率は望めないことである。そこで、一般化した階層的な並列処理を提案している。これはプロセッサを幾つかのグループに分け、1 グループを 1 つの並列計算機に見立てて、その中でフラグメントの計算を並列に実行する方法である。グループ数およびグループ内のプロセッサ数の割り振りを、フラグメントの数やサイズをもとに計算時間が最短となるよう適化している。
幾つかの QCLO の応用については本討論会にて別途報告するので、併せてご参考頂きたい。
本研究は文部科学省 IT プログラム「戦略的基盤ソフトウェアの開発」において実施された。
【参考文献】[1] Sato, F. et al,, Int. J. Quant. Chem. 63 (1997) 245-256; Sato, F. et al, Chem. Phys. Lett. 341 (2001) 645-651.
[2] Kashiwagi, H. et al, Mol. Phys. 101 (2003) 81-86.

1Dp07

大規模タンパク質全電子計算を目指した
地球シミュレータ版 ProteinDF の開発
吉廣保†、稲葉亨‡*、上野哲哉†、佐藤文俊†、柏木浩†‡

東京大学生産技術研究所†

アドバンスソフト株式会社‡株式会社日立製作所*

１．はじめに
生体内で重要な働きをするタンパク質は多くの原子からなる生体高分子である。特に大型のタンパク質においては、基本骨格を作る 20 種類のアミノ酸以外にも、ヘムやクロロフィルに代表されるような色素分子が複雑に入り込んで、機能性に富んだタンパク質固有の反応を発現している。タンパク質の機能を理論的に電子レベルで解明するためには、タンパク質全体を分子として扱うことができ、かつ多くの色素に含まれる金属もアミノ酸残基と同程度の精度で求めることができる解析手法が望ましい。
そのような観点から、発表者らはこれまで独自にガウス型基底関数を用いた密度汎関数法に基づく大規模分子軌道法プログラム ProteinDF を作成・発展させてきた[1]。密度汎関数法は
Kohn‑Sham 方程式に基づく分子軌道法の一種で、最大の利点は ab initio HF 法と同程度の計算量で、電子相関効果を取り込んだ計算ができることである。そのため、大型のタンパク質の計算に威力を発揮するであろう。また、開発したプログラムは、カノニカル分子軌道計算においても、計算時間が軌道数の約 2.4 乗で計算できるという良好な結果を得ている[2]。
２．オブジェクト指向技術による並列化と WS クラスタによる並列分散処理
ProteinDF プログラムは大規模な行列等のデータを安全に取り扱い、かつプログラムを多人数で共同開発することを念頭に置き、オブジェクト指向言語 C++でコーディングされている[1]。
発表者らは、オブジェクト指向の考え方を並列化にも導入して、プログラム全体をアルゴリ
ズム、計算、通信といった各機能が互いに独立した構造となるよう設計したうえで、シナリオ（アルゴリズム）部、(並列) 計算部、仮想マシン部、通信部に完全に階層化し、並列化した(図1)。これらの階層は、機能的に全て独立しているので、計算機アーキテクチャや通信用ライブ
ラリの多様性に容易に対応できる。また、各計算ルーチンが十分に独立しているため、他の計算ルーチンのデータ構造などを考えることなく個
別に修正・改良を加えることができる。
Master Process
Scenario Object
こ の 並 列 版 ProteinDF を 用 い て 、 HP ALPHA
21264 667MHz ワ ー ク ス テ ー シ ョ ン （ WS ） と 同
Parallel Computational Object
21164A 533MHz WS 8 台、500MHz WS 6 台の計 15 台
を 100Base‑TX の Ethernet で接続した WS クラスタ
Virtual Machine Object 1
Virtual Machine Object n
並列システムを使用し、馬心臓シトクロム c（104Communication
Communication…
ObjectObject残基+ヘム c、1738 原子、9600 軌道、17578 補助基
PVM
PVM底関数）の全電子計算に成功した実績がある。しかし、このシステムでは繰り返し計算の 1 回にかPVM
PVMCommunication
Communication
かる実時間は 77362 秒(21.5 時間)で、59 回の繰り…
ObjectObject
返し計算を要したため[2]、それほど手軽に実行でVirtual Machine Object
Virtual Machine Object
きるわけではない。また、さらに大きなタンパク
ActualActual質の計算を行うには、この程度の計算機資源ではComputational Object
Computational Object
不十分である。
Slave Process 1

３．地球シミュレータへの移植とチューニング

Slave Process n

図 1：ProteinDF の並列化構造

そこで、本研究では地球シミュレータのような超並列大型計算機において ProteinDF が効率良く実行できるプログラムを開発している。本研究により、100 残基から数百残基のタンパク質

の全電子計算が、現実的な時間で実行できる基盤が整備できるものと考えられる。
地球シミュレータは、640 台の計算ノードを、640×640 の単段クロスバネットワークで結合させた分散メモリ型並列計算機である。各計算ノードは、ピーク性能 8Gflops のベクトル型計
算プロセッサ 8 台が主記憶装置 16GB を共有する共有メモリ型並列計算機となっている。全体では計算プロセッサが 5120 台でピーク性能は 40Tflops、主記憶容量は 10TB となる計算機である。
このような計算機資源で ProteinDF が効率よく実行できるかを図る目安として、地球シミュ
レータに比べ、CPU の性能に対して遥かに低速なネットワーク性能を持つ WS クラスタ並列シス
テムにおいて、どの程度の効率で ProteinDF が並列化されているかを調べた。計算には 31 残基
ペプチド鎖断片（2867 軌道、5273 補助基底関数）の計算を使用した。HP Alpha 21164A 533MHz
WS 8 台を 100Base‑TX の Ethernet で接続した 1,4,8 台のシステムを使用して、最初の SCF 計算にかかる実時間を測定し、1 台で計算した時間を基に使用台数による並列化効率を求めた。その結果、並列化効率は 4 台で 91.0%、8 台で 78.1%であった[3]。解析の結果、並列化効率を下げている最も大きな要因は行列データの転送時間であり、計算処理は 95%以上の効率で均等化されていることが明らかとなった。行列の転送時間よりも計算処理時間のほうが計算規模に対する依存性が大きいため、さらに大きな系ではより高い効率が見込まれる。これらの結果により、地球シミュレータにおいて高速なシミュレーションが達成できるものと考えられる。
地球シミュレータでの並列処理は、OpenMP や Message Passing Interface (MPI)を使用して並列化を行う。しかし、OpenMP も MPI も SPMD 方式の並列処理に優れた処理環境やライブラリである。並列版 ProteinDF は巨大計算のためにメモリをより柔軟に利用できるよう、マスタ・ス
レーブ方式を採用している。マスタ・スレーブ方式でプログラムをコントロールするためには
動的プロセス制御を行える通信ライブラリが適しており、これまで Parallel Virtual Machine(PVM)を採用してきた。地球シミュレータにおいても、数万次元の行列を扱うには、計算ノード
当たり 16GB のメモリを効率よく使用する必要があるので、本研究では、ProteinDF の地球シミ
ュレータへ移植にも、あえてマスタ・スレーブ方式を採用し、動的プロセス生成のできる MPI‑2
通信ライブラリを用いた。
通信クラスの MPI‑2 通信ライブラリへの対応の際、通信クラス以外の部分においても一部 MPI
通信でのメッセージサイズを明示的に指示するための変更が必要となったが、図 1 のプログラ
ム構造のお陰で、通信クラス以外にはほとんど手を加えること無く、PVM 版通信クラスと MPI‑2
版通信クラスのどちらでも動作するプログラムを作成することが出来た。また、本研究における、現行の地球シミュレータの利用形態は、計算ノード内・計算ノード間を特に区別せずすべての計算プロセッサを対等に扱った、first‑come/first‑serv タイプのスレーブ計算処理スケジューリングによる並列処理法を採用している。
地球シミュレータはベクトル型計算プロセッサの並列マシンである。ProteinDF は分子積分や行列演算などが律速ルーチンである。各ルーチンは効率よく並列化されたので、続いてスレー
ブプロセッサで実行される実計算ルーチンのベクトル化を行っている。現在のところ、分子積
分ルーチンはベクトル化率約 92%で、行列計算は専用ライブラリをスレーブプロセスで用いることで性能が得られている。
発表では、さらにチューニングを施した上、シトクロム c の全電子計算について WS クラスタと地球シミュレータでの計算結果をもとに、並列化効率などについて報告する予定である。
謝辞
本研究は文部科学省 IT プログラム「戦略的基盤ソフトウェアの開発」において実施された。
さらに、地球シミュレータセンターに感謝します。
参考文献[1] F. Sato, Y. Shigemitsu, I. Okazaki, S. Yahiro, M. Fukue, S. Kozuru, H. Kashiwagi, Int. J. Quantum Chem.
63, 245 (1997).
[2] F. Sato, T. Yoshihiro, M. Era, H. Kashiwagi, Chem. Phys. Lett. 341, 645 (2001).
[3] T. Yoshihiro, F. Sato, H. Kashiwagi, Chem. Phys. Lett. 346, 313 (2001).

1Dp08

蛋白・核酸理論研究のための非経験的量子化学計算ツールの開発(豊橋技科 a 科技団ｂ)○ 関野 秀男 a ・杉木真一郎 a,b・仙石康雄 a・井上英計 a・澤田寛明 a・栗田 典之 a

巨大生体高分子などナノスケールの非経験的量子計算では全系の完全に変分的記述は現段階ではかなり困難であり、また生体活性の理論解明という観点からはあまりに些細な artifact 解決にほとんどの計算労力を費やす結果になりかねない。
我々のグループではこうしたシステムの理論解明と技術的制約とのトレードオフをはじめから念頭に入れ、プログラム開発を進めている。蛋白や核酸などの生体高分子はその構造の柔軟性と、いたるところに存在する電気的分極により生体機能の発現を誘発する。分極をとおした構成フラグメント間の相互作用は弱いが多くの箇所で同時に働き有機的に作用する。こうした相互作用を量子論的に表現するのには多体効果である電子相関が重要となってくる。然し伝統的な量子化学の方法論による
電子相関効果の見積りは多大な計算労力を要し、現在のところこうした巨大系への応用をめざしたプログラムの基礎アルゴリズムとして採用するのには多くの乗り越えるべき困難を抱えている。我々は密度汎関数法(DFT)をフラグメント分子軌道法(FMO)の枠組みに採り入れることによりプログラム開発をおこなっている。多くの化学結合の理論解明に有力である汎関数が分子間力の記述に無力であることから、汎関数の選択や改良も同時に行っている。DFT FMO 法は Hartree Fock (HF) FMO
法のほぼ延長上にあるが HF FMO 法では起こらなかった、フラグメント電荷の収束困難という新たな問題に遭遇した。

現在そうした問題を解決し、並列化により現

実系への応用が可能なプログラムが完成しつつある [4]。FMO 法による全エネルギ
ーが近似なしの計算に比べかなり定量的に満足の行くものである事は本方法の提唱者である北浦[1,2]や他の日本人研究者等[3]により示されてきたが、我々は主に生体分子の活性やナノスケール分子の物性算定に関する定性的理論枠組みとしてどれほ
ど有用であるがと言う観点で我々の理論とそれに基づいて構築されたプログラムの数値的検証を行ったので報告する。
FMO の全エネルギー算出に用いる密度行列を用い一次の物性を算定

P = ∑ PIJ − ( N − 2)∑ PI
I >J

PI = Tr PI ρ I

I

PIJ = Tr PIJ ρ IJ

すると近似なしの 算定値をかなり良く再現することから、ゼロ次の FMO とそのペア補正がかなり良い近似であることがわかるが[5]、更に高次の物性もかなり良く再現できることが分かった。 計算の些細や FMO が良い近似とならない応用例などについては当日発表する。

[1] K.Kitaura, T.Sawai, T.Asada,T.Nakano,M.Uebayashi, Chem. Phys. Lett. 312, (1999) 319.
[2] K.Kitaura, E.Ikeo, T.Asada, T. Nakano, M. Uebayashi, Chem. Phys. Lett. 313, (1999) 701.
[3] T. Nakano, T.Kaminuma, T. Sato, Y. Akiyama, M. Uebayashi, K. Kitaura, Chem. Phys.
Lett. 318, (2000) 614.
[4] S. Sugiki, Y. Sengoku, N. Kurita and H. Sekino, submitted.
[5] H . Sekino, Y. Sengoku, S. Sugiki and N. Kurita, to appear in Chem. Phys. Lett..

本研究は、JST-ACT (研究開発課題「DNA のナノ領域ダイナミクスの第一原理的解析」) の援助を受けて行われた。