Gaussian教程–ECD计算
摘要：在天然产物的结构确证中经常需要用电子圆二色性（electronic circular dichroism，ECD）来确定化合物绝对构型，通过TD-DFT理论计算ECD图谱与实验ECD图谱比较可以辅助立体化学结构确证。本教程以Gaussian软件为基础，描述了TD-DFT理论计算ECD的详细过程。
作者：肖高铿
电邮: gkxiao@molcalx.com
一. ECD计算简介
确定化合物的绝对构型（absolute configuration，AC)在天然产物的结构确证里经常是个挑战性的工作。当一个手性化合物具有合适生色团，电子圆二色性（electronic circular dichroism，ECD）就或许可以成为结构确诊的强力工具。近几年，越来越多的人用TDDFT（time-dependent density functional theory）理论来进行ECD计算辅助结构确证。也越来越多的综述对ECD计算进行总结，比如Nugroho A. E.与H. Morita(2014)就详细地总结、评述了ECD计算的方法与流程；另外《国际药学研究杂志》2015年12月第42卷第6期为一专刊《解释天然产物立体结构的技术与方法》，本期就有很多文章对ECD计算进行结构确证进行了细致的描述。
本教程按照Nugroho, A. E.与H. Morita(2014)的计算流程，详细描述如何进行计算预测理论ECD光谱。
二. ECD计算流程
ECD计算流程如图1所示。假设一个化合物是柔性的，那么它在溶液中就存在多个构象，实验测得的UV、ECD是这个化合物多个构象在时间与空间的平均，因此需要对这个化合物进行构象搜索获得低能构象集，再分别对每个构象进行UV与ECD计算，最后采用玻尔兹曼平均获得最终的图谱。
图1. ECD计算流程图
第一步. 构象搜索
流程中的构象搜索可以用CONFLEX或Torch来实现。CONFLEX是一款用户非常友好的穷尽构象搜索软件，可以非常方便的完成构象搜索。它采用MMFF94S力场，报告每个构象占据构象空间的百分比，可以依此比例挑出合理的低能构象集用于进一步的结构优化；而Torch也是一款用户非常友好的构象搜索软件，在构象搜索前可以设定能量窗、并去除重复构象以控制构象多样性（通过设定RMS阈值来控制）。
CONFLEX构象搜索教程：http://blog.molcalx.com.cn/2016/04/12/conflex-conformation-search.html
Torch构象搜索教程：http://blog.molcalx.com.cn/2016/06/04/torch-tutorial-conformation-hunting.html
第二步. 采用量子化学计算对第一步获得的每个构象进行再次优化
第一步采用分子力学获得了低能构象集，排除了不可能的构象；在第一步的基础上，需要对每个构象分别采用DFT进行再次优化、并将溶剂考虑进来以获得真正的低能构象。本步计算采用Gaussian 09/16软件来进行。
第三步. 采用TDDFT计算激发态
在第二步优化的低能构象基础上，采用相同的基组与方法进行TDDFT计算。本步计算采用Gaussian 09/16软件来进行。
第四步. 数据处理与图谱生成
SpecDis是Bruhn T(2013)等人专为生成VCD、ECD、UV等图谱而设计的软件。SpecDis可以直接读入第三步的Gaussian计算结果并画出UV与ECD图谱；还可以将几个图谱进行玻尔兹曼平均获得最终图谱；还可以将计算的图谱与实验图谱进行定性、定量地比较以辅助确定绝对构型。
GaussView是Gaussian配套的可视化软件，可以用来准备计算文件、作业管理与结果分析。对于ECD计算，可以用来生成UV与ECD图谱。GaussView 6支持对多个构象的图谱进行玻尔兹曼加权模拟新的图谱，但是GaussView 5.08或以前版本并不支持该功能。
总的来说，ECD计算可以分为四个步骤来完成，分别用到不同的软件，总结如图3所示。
图3. ECD计算步骤与用到的软件
三. 操作
因为CONFLEX已经另有教程且SpecDic的使用极其简单，所以本教程主要提供了第二、三步的操作。假设CONFLEX计算后，每个构象被保存为一个sdf格式文件，其中第1、2、3...构象分别保存为CONF_1.sdf, CONF_2.sdf, CONF_3.sdf...。我们现在以其中的第一个构象CONF_1.sdf为例说明如何进行上述的量化计算再优化与激发态计算。
1. 采用量子化学计算对第一步获得的每个构象进行再次优化
(1) 用GaussView读入CONF_1.sdf
File | Open
在Open File对话框里设定文件类型为MDL FILES, 选取CONF_1.sdf, 点击打开按钮。
注：GaussView支持Mol2，Mol，SDF，PDB等常见的文件格式。
(2) 生成作业文件
Calculate | Gaussian Calculation Set up
弹出Gaussian Calculation Setup对话框，设定Job Type为Optimization；
点击Edit, 弹出Gaussian Calculation对话框，点击Save，将作业保存为CONF_1_min.gjf；此时会弹出写字版编辑器，编辑以#开头的Route部分，确保其内容为：
1 | #p opt apfd/6-311+g(2d,p) SCRF=(IEFPCM,Solvent=Methanol) integral(ultrafinegrid) |
点击保存，关闭写字板。
(3) 提交作业
关闭写字板后会弹出Run Gaussian对话框，点击OK开始进行计算。也可以取消作业，之后再用Gaussian打开上述的作业文件CONF_1_min.gjf进行计算。
Linux用户在shell里键入命令：
1 | g09 CONF_1_min.gjf |
计算完毕获得CONF_1_min.chk, CONF_1_min.log两个文件。
(4) 小结
总的来说，在这个练习中，我们用GaussView读入了CONFLEX搜索获得的结构，并用GaussView生成了结构优化的作业文件，并进一步提交了该作业。本次的结构优化在apfd/6-311+g(2d,p)理论水平进行，并用PCM溶剂模型考虑了甲醇的溶剂效应。
注意：建议追加个频率计算以确认能量最小化计算获得的结构不含虚频，Route部分：
1 | #p opt freq apfd/6-311+g(2d,p) SCRF=(IEFPCM,Solvent=Methanol) integral(ultrafinegrid) |
2. 激发态计算(UV与ECD)
用GaussView打开上一步生成的log文件CONF_1_min.log，用GaussView读入该文件，用类似上步优化的作业文件过程，新建一个作业存为CONF_1_ecd.gjf，Route部分内容如下：
1 | #p td=(nstates=30) apfd/6-311+g(2d,p) SCRF=(IEFPCM,Solvent=Methanol) integral(ultrafinegrid) |
在GV里可直接提交任务，而Linux用户可在shell里键入命令：
1 | g09 CONF_1_ecd.gjf |
计算获得CONF_1_ecd.chk与CONF_1_ecd.log两个文件。
注：
- 1）上面推荐用的apfd/6-311+g(2d,p)是Gaussian 09 D.01之后推荐的标准方法，与B3LYP相比，APFD增加了色散作用项，但是计算速度与B3LYP相当。
- 2）nstates=30，视实验图谱与计算图谱的差异而增减，如果实验图谱比计算图谱多出峰，请增加。
3. 图谱处理
TD计算给出State的波长与振子强度、转子强度，将波长~振子强度画图，得紫外图谱；波长~转子强度画图，获得得ECD图谱。如何用Excel作图，可具体参考Gaussian 网站教程。更加简便的处理方法是用GaussView或第三方软件SpecDis等来进行处理。
用GaussView或SpecDis读入CONF_1_ecd.log文件(见图4)，就可以生成相应的UV与ECD图谱（见图5），可以将数据导出到Excel里与实验值进行比较（图6），这里不再详细进行叙述。下面的视频演示了如何用GaussView生成图谱的过程。在这个视频中，两个图谱根据能量计算的玻尔兹曼权重混合出新的图谱。当然，权重还可以根据频率计算的构象自由能来计算，具体方法参见NMR计算教程。也演示了如何调整半峰宽以调整峰形，以便使得计算图谱与实验图谱一致。
需要提醒的是，计算值与实验值的峰位有一定的系统误差，一定要使用UV进行shift值（偏移量）校正。计算的UV与实验UV峰位进行比较，确定偏移量。如果不使用UV提供的校正，SpecDis经常会给出一个很大Shift值的ECD图谱，这种情况下科顿效应是不能得到正确归属的。出现这个现象的原因是：SpecDis用曲线下面积差（越小越好）作为打分函数，导致-50nM的偏移比之-3nM的偏移打分仅仅好那么一点点SpecDis就将-50nM的偏移作为最佳结果推荐给用户。过高的UV shift值意味着过拟合，此情况下得到的结论必然不靠谱，阅读文献时要非常小心，很可能文献定的结构就是错误的。
图4. GaussView可以读入激发态计算结果，绘制UV与ECD图谱
图5. GaussView可以读入激发态计算结果，绘制UV与ECD图谱
图6. GaussView可以将UV,ECD等图谱以数据点格式导出给第三方软件使用，比如Excel：在图谱上点击右键>Save Data
图7. Excel模拟一个化合物8个构象的图谱及其玻尔兹曼平均图谱
四. 用Excel处理ECD图谱
从GaussView里导出计算的UV或ECD图谱的数据点之后（见图6），可以读入Excel与实验图谱进行比较。最常问的问题是如何平移计算的图谱以便与实验图谱进行比较。下面分几个小点来讨论一下。
- 确定计算图谱要平移的大小，即偏移量
- 曲线平移预备知识
- 平移与缩放曲线
- 如何制作散点图
推荐根据实验值来设定偏移量：根据实验测得的UV吸收峰与对应计算图谱的UV吸收峰进行比较，得到一个差值，该差值就是偏移量。将ECD谱按照该偏移量进行平移校正。往那边移动，最笨的办法是自己试试就知道了。
对于y=f(x)这样函数的曲线平移，在高中的教科书有详细论述。我百度了一个《高中数学平面曲线平移伸缩变换的技巧》[6]，详细地介绍了进行平面内曲线的平移与伸缩变换的数学原理，请务必掌握，因为等下需要用该知识来指导下面的excel操作。
无论是计算的UV或ECD图谱，还是实验的UV或ECD图谱，都至少包含两列：一是作为横轴（x）的波长，二是作为纵轴（y）的Epsilon。根据曲线偏移预备的知识，平移一个曲线是保持纵轴不变，而对横轴的数值进行调整。比如，向左平移5nM，就是将波长列减去5nM，再绘制散点图。
同样的道理，缩放曲线（调整曲线高度）是保留横轴（x轴）不变而对纵轴（y轴）进行变化。那么只要新增一列，该列是对原有y轴的数值的调整，比如乘以一个系数，然后再绘制散点图。
具体的操作，请参考百度经验《Excel 技巧：如何上下左右平移图形》。
请参考百度经验《excel散点图怎么做》
相关主题
- Gaussian教程 | 计算分子的拉曼光学活性(ROA)
- Gaussian教程 | 计算分子的比旋光度
- Gaussian教程 | VCD计算
- Gaussian教程 | NMR图谱与耦合常数的计算
- ECD图谱预测与生成Excel报告
五. 文献
[1] Nugroho, A. E. and H. Morita (2014). Circular dichroism calculation for natural products. J Nat Med 68(1): 1-10. DOI：10.1007/s11418-013-0768-x
[2]Bruhn, T.; Schaumlöffel, A.; Hemberger, Y.; Bringmann, G.; SpecDis Version 1.64 ed.; University of Würzburg: Würzburg, Germany, 2013.
[3]Bruhn, T.; Schaumlöffel, A.; Hemberger, Y.; Bringmann, G. SpecDis: Quantifying the Comparison of Calculated and Experimental Electronic Circular Dichroism Spectra, Chirality, 2013, 25, 243–249.
[4]Gaussian Inc.(2016,February 25) Ploting UV/VIS SPECTRA FROM OSCILLATOR AND DIPOLE STRENGTH. Retrieved from http://www.gaussian.com/g_whitepap/tn_uvvisplot.htm
[5] 陈国栋，肖高铿，姚新生，高昊. 电子圆二色谱计算方法在天然产物绝对构型确定中的应用. 国际药学研究杂志. 2015, 42(6):738-743
[6] 蔡正伟. 平面内曲线平移伸缩变换的技巧.百度文库：高中数学平面曲线平移伸缩变换的技巧. 2019年12月3日访问：https://wenku.baidu.com/view/21ac76c5541810a6f524ccbff121dd36a32dc480.html
六. 联系我们
电邮：info@molcalx.com
电话：020-38261356
单位：广州市墨灵格信息科技有限公司
网站：http://www.molcalx.com.cn/gaussian
CONFLEX软件试用申请与下载：http://www.conflex.net
咨询ECD全套方案：软硬件配置推荐
原创文章，作者：小墨，如若转载，请注明出处：《Gaussian教程–ECD计算》http://blog.molcalx.com.cn/2016/05/23/gaussian-ecd-tutorial.html
