幾種生成有機分子GROMACS拓撲文件的工具
幾種生成有機分子GROMACS拓撲文件的工具
文/Sobereva @北京科音
First release：2014-Dec-9 Last update: 2022-Mar-15
本文把各種可以產生GROMACS拓撲文件的工具進行匯總。具體細節和實際使用例子筆者在北京科音分子動力學與GROMACS培訓班里會詳細講（http://www.keinsci.com/workshop/KGMX_content.html）。
值得一提的是，這些程序有很多不能產生原子電荷，或者產生的原子電荷質量較差。RESP電荷是最適合分子動力學模擬使用的電荷之一，AMBER、GAFF、GLYCAM等力場將之作為御用的原子電荷，將RESP電荷與UFF力場結合使用也是得當的。Multiwfn是計算RESP電荷最簡單、最好、最靈活而且完全免費的工具，原理、用法和實例見《RESP擬合靜電勢電荷的原理以及在Multiwfn中的計算》（http://www.shanxitv.org/441）、《計算RESP原子電荷的超級懶人腳本（一行命令就算出結果）》（http://www.shanxitv.org/476）、《RESP2原子電荷的思想以及在Multiwfn中的計算》（http://www.shanxitv.org/531）。另外，OPLS-AA力場很適合結合1.2*CM5電荷使用，用Multiwfn的主功能7的子功能16直接就能計算出CM5電荷，將之手動乘以1.2即是1.2*CM5電荷，可以用Gaussian、ORCA等程序產生的fch/wfn/wfx/molden等格式作為Multiwfn的輸入文件，詳見《詳談Multiwfn支持的輸入文件類型、產生方法以及相互轉換》（http://www.shanxitv.org/379）。如果你完全不會用Gaussian的話，可以直接用這個傻瓜式腳本，一行命令就能算出來：《計算適用于OPLS-AA力場做模擬的1.2*CM5原子電荷的懶人腳本》（http://www.shanxitv.org/585）。另有專門結合ORCA用的傻瓜式腳本，見《ORCA結合Multiwfn計算RESP、RESP2和1.2*CM5原子電荷的懶人腳本》（http://www.shanxitv.org/637）。
? Sobtop（http://www.shanxitv.org/soft/Sobtop）
這是我開發的GROMACS拓撲文件產生工具，主要產生GAFF、AMBER力場的拓撲文件，但由于其力場庫可以自行非常方便地修改和擴充，因此Sobtop本質上是完全普適、通用的。Sobtop可謂是最理想、最靈活、最易用的產生GROMACS拓撲文件的工具。此程序用起來超級簡單，什么額外的程序以及特殊的運行環境都不需要裝，解壓即用。Sobtop使用極其方便，照著屏幕上的提示敲幾下鍵盤，itp、top和gro文件就產生了，另外也可以要求產生rtp文件。Sobtop的主頁有非常詳細的產生各類體系拓撲文件的例子，并給出了詳細的相關要點的說明。從例子中你會體會到Sobtop的設計特別注重兼顧便利和靈活，初級用戶會體會到它極其便利，而高級用戶則會體會到通過Sobtop構建復雜體系拓撲文件特別靈活好用。
Sobtop可以產生任意體系的拓撲文件，有機和無機小分子、過渡金屬配合物、聚合物、共價團簇、晶體、二維材料等等全都可以產生，孤立體系和周期性體系都支持。Sobtop可以用于包含任意元素的體系，那些GAFF/AMBER不支持的元素可以自動讓Sobtop用UFF的，或者自己定義額外原子類型。對于GAFF/AMBER力場里缺失的成鍵相關參數，可以直接讓Sobtop基于量子化學程序產生的Hessian矩陣通過不同方法自動算出來，也可以自己把其它方式得到的（如自己單獨擬合、文獻中找的）添到力場庫文件里讓Sobtop直接用。Sobtop運行效率特別高，甚至對于幾千原子體系的拓撲文件都可以搞。
Sobtop極度靈活，提供許多不同的指認原子類型的方式，可以自動指認GAFF/AMBER原子類型，也可以手動指認，也可以自定義判斷規則讓Sobtop結合實際化學環境和成鍵關系判斷，等等。Sobtop還提供了豐富的指認力場成鍵相關參數的方式，比如可以都用力場庫文件里的，都用直接基于Hessian算出來的參數，或者一部分作為剛性（里面的參數都基于Hessian算出來的）而其余部分作為柔性（參數用力場庫里的，可以考慮二面角的可旋轉性），等等。Sobtop還可以調用Multiwfn或OpenBabel自動指認EEM、Gasteiger、MMFF94原子電荷，還可以載入Multiwfn計算RESP/RESP2等電荷產生的chg文件來獲得原子電荷。Sobtop的設計在各方面都特別貼心，比如在產生的拓撲文件里非常詳細地注釋各個力場項的參數來源、是否是缺失的，在自己手動改和補參數時很方便。
有了Sobtop，下面介紹的acpype就沒有任何使用價值了。acpype只適合產生GAFF能描述的那些有機和極少數無機體系，有特殊成鍵方式、有GAFF不支持的元素、周期性體系都沒法搞；對很大體系耗時非常高；容易出現莫名其妙又不好解決的報錯；還得有Python運行環境、安裝臃腫的AmberTools...
? acpype：這是一個Python腳本，可以在http://svn.code.sf.net/p/ccpn/code/branches/stable/ccpn/python/acpype/下載。使用前必須在機子里安裝AmberTools（免費），acpype會調用其中的Antechamber先產生Amber格式的拓撲文件然后再轉成GROMACS的。acpype用法很簡單，要處理xxx.mol2就執行./acpype.py -i xxx.mol2，算完后會新產生一個xxx目錄，里頭有_GMX后綴的.gro、.itp、.top，直接在GROMACS里用即可。默認情況下，產生的拓撲文件是基于GAFF力場的，另外也會輸出_OPLS后綴的基于OPLS力場的文件，但屬于實驗性質不建議用。.mol2文件用常用的GaussView就可以產生（但必須確保在gview里看到的分子結構中沒有詭異的成鍵方式），也可以通過OpenBabel或Antechamber將其它格式轉成.mol2。默認情況下acpype分配的原子電荷是Antechamber產生的AM1-BCC，雖然能用，但明顯不如RESP/RESP2電荷理想。建議大家按前述做法用Multiwfn計算出RESP或RESP2電荷，自行寫入分子拓撲信息的[atoms]的原子電荷那一列。
如果你懶得為了用acpype而裝臃腫龐大的AmberTools，可以用在線版http://bio2byte.be/acpype/，不過可能排隊要排很久。所以如果你要快速處理較多小分子，還是建議用離線版acpype。
對于稍大的體系，用獨立的acpype時強烈建議加上-c user選項以避免acpype自動算AM1-BCC原子電荷，否則在處理期間acpype所調用的Antechamber會先調用SQM程序用半經驗方法進行優化然后再算這個電荷。然而SQM不僅慢，做優化還容易不收斂，導致半天也無法成功產生拓撲文件。更何況最后也是要替換為Multiwfn算的RESP/RESP2電荷，故計算AM1-BCC電荷沒實際意義。同理，用在線版acpype也是建議把charge method設為user。
? PRODRG2（http://davapc1.bioch.dundee.ac.uk/cgi-bin/prodrg）：歷史非常悠久很有名的生成Gromacs拓撲文件的在線工具。只支持GROMOS87/96力場，生成gromacs和其它一些程序的拓撲文件。原子電荷是根據基團指認的，如果基團識別不對那么電荷也不可靠。有在線版和離線版。可以選擇自動優化分子結構。此工具在J. Chem. Inf. Model., 50, 2221 (2010)被曝光往往不能正確指認charge group，導致原子電荷分配也很不合理。而且又由于有了更好的ATB，因此如今強烈不推薦使用。
? Automated Topology Builder（ATB，http://atb.uq.edu.au）：生成GROMACS拓撲文件的在線工具，可以生成基于ATB開發者自行修改的GROMOS96 G54A7力場（擴充了原子類型）的GROMACS或GROMOS程序的拓撲文件。比PRODRG2更先進更可靠，解決了PRODRG2沒法自動確定質子化態和charge group指認不準的問題以保證原子電荷可靠，并且能利用對稱性保證電荷等價，另外沒有PRODRG2那樣對于每日提交的數目有限制。ATB雖然設計得不錯，可以作為產生小分子GROMOS力場的GROMACS拓撲文件的首選，但它生成的原子電荷、確定的參數的可靠性也只是一般，不能保證總是很理想。如果要做很嚴謹的模擬研究，還是建議自行計算RESP電荷并且手動檢查ATB給出的拓撲文件的合理性，有必要時適當調節。ATB網站上也有個分子庫，包含巨量事先搞好的小分子的參數和拓撲文件，其中有的是別人之前提交過的分子，有的是經過專人手工處理過的分子，顯然后者可靠程度更高。對于比較常見的分子，提交到ATB之前建議先搜索一下分子庫，如果已經有的話就直接用。
? MKTOP：是個Perl腳本，可以在 注：以前MKTOP是有在線服務器的（http://www.aribeiro.net.br/mktop/），但是如今已經失效。
? TPPmktop（http://erg.biophys.msu.ru/tpp/）：在線工具，提供pdb文件，能產生OPLS-AA力場參數的GROMACS拓撲文件。速度比較快，但得到的拓撲文件里有時候會缺參數，或者出現額外的原子類型，需要再手工處理（需要引入額外的力場文件，但是筆者發現在網上下載不到）。此服務器有時候有其它任務在跑，此時無法提交任務，只能等過一陣服務器沒任務時再提交。
? LigParGen（http://zarbi.chem.yale.edu/ligpargen/）：生成GROMACS, NAMD, CHARMM, LAMMPS等程序OPLS-AA力場的拓撲文件的工具。畢竟這是OPLS力場開發者自己搞的，因此原則上比MKTOP、TPPmktop更靠譜。可以通過SMILES字符串、.mol、.pdb文件進行輸入（用pdb格式時老是報錯，我建議用.mol），原子上限200個，可以順帶著讓服務器對分子在OPLS-AA力場下做幾何優化。分配的原子電荷是1.14*CM1A或1.14*CM1A-LBCC，前者是把CM1A電荷數值乘上1.14得到的，后者是在1.14*CM1A基礎上再引入LBCC校正得到的，在J. Phys. Chem. B, 121, 3864 (2017)中已證明這兩種原子電荷結合OPLS-AA模擬有機體系凝聚相可以得到不錯結果。此服務器還可以輸出PQR文件。
本文開頭所述的Multiwfn可以計算的1.2*CM5電荷比LigParGen直接給出的電荷在多數情況下更適合做動力學模擬，特別是對于密度、蒸發焓的精度方面而言，這點在J. Phys. Chem. B, 121, 3864 (2017)的測試中充分體現了。而且對于一些分子，LigParGen給出的原子電荷加和不精確為0，而是比如0.0001，這導致此類分子數目較多時體系總電荷對整數有不可忽視的偏離，對模擬造成不良影響。因此建議大家將LigParGen給出的itp文件里的原子電荷部分替換為Multiwfn得到的1.2*CM5電荷。
? OBGMX（http://software-lisc.fbk.eu/obgmx/）：生成UFF力場的gromacs拓撲文件的在線工具。由于UFF幾乎涵蓋整個周期表，因此不光有機分子，也可以使得Gromacs能夠處理無機物、周期性體系，比如MOF。由于UFF的力場形式和Gromacs所支持的不完全兼容，此程序給出的參數實際上是對原UFF力場的近似。個別原子類型識別可能有誤。電荷必須自己提供。
2020-Dec-26注：此在線服務器目前已無法使用，但可以下載離線程序使用，見http://bbs.keinsci.com/thread-21015-1-1.html。
? SwissParam（http://swissparam.ch）：輸入有機小分子mol2文件，生成用于CHARMM/NAMD和GROMACS模擬的拓撲、參數文件。力場參數基于MMFF，但只保留諧振項部分，因此只是MMFF的近似。原子電荷通過MMFF方法獲得。范德華參數采用CHARMM22中最接近的原子類型。這樣的參數比較粗糙，有優化的余地。
? CGenFF（https://cgenff.paramchem.org）：先注冊，上傳有機小分子的mol2文件，即可生成基于CGenFF力場的CHARMM的拓撲文件。如果mol2是gview建的，一定要事先把里面的Ar替換成ar。上傳文件的時候不要選Guess bond orders from connectivity和Include parameters that are already in CGenFF。如果文件處理正常，會立刻產生str文件，點擊其鏈接之后把里面內容拷到比如Actos.str里面。進入網頁里的More Info & Tools - Utilities，點擊GROMACS conversion program，可下載把str文件轉換為GROMACS格式的Python腳本cgenff_charmm2gmx.py。對于CentOS 7.2，應運行yum install numpy和yum install python-networkx把這個腳本所需的包裝上。去http://mackerell.umaryland.edu/charmm_ff.shtml#gromacs里面下載用于GROMACS的CHARMM36力場文件，解壓得到比如charmm36-nov2016.ff文件夾，將它和Actos.str、cgenff_charmm2gmx.py都放到當前目錄，另外也把此文件夾拷到gromacs的top目錄下。假設Actos.str里RESI后面的詞是Molecu，最初上傳的是Actos.mol2，則運行./cgenff_charmm2gmx.py Molecu Actos.mol2 Actos.str charmm36-nov2016.ff。在當前目錄會產生molecu.top、molecu.prm、molecu.itp、molecu.pdb（從mol2轉換過來的），適當調整itp和top里的分子名，調整itp里的殘基名和結構文件相對應后，即可用于模擬。str文件里的力場參數有penalty指標，數值越大說明此參數可靠性越低，詳見網頁里的說明。
? ztop：計算化學公社論壇（http://bbs.keinsci.com）上鐘成開發的拓撲文件產生工具，請在論壇首頁搜索框里搜索ztop看他發的相關帖子了解此程序的特點和使用。
另外說兩個有關的程序
? YASARA AutoSMILES Server（http://www.yasara.org/autosmilesserver.htm）是在線版程序，離線的話得買YASARA（建模+可視化+動力學模擬工具）。可以計算AM1-BCC電荷、指認amber94/96/99力場參數。但是產生的文件只能用YASARA View打開，也就是說，拓撲文件相當于YASARA專用的。
? RED（RESP ESP charge Derive，http://q4md-forcefieldtools.org/REDServer-Development/）：專門生成RESP電荷的，也能搞力場參數。頁面做得不是一般的糟糕，結構十分混亂，令我摸不到頭緒，因此筆者也沒怎么仔細研究。由于Multiwfn已經完美支持RESP電荷計算了，RED已沒有必要使用了。
這里把各種程序做一下總結，是筆者講授的北京科音分子動力學與GROMACS培訓班中的一頁幻燈片
