2019年4月29日-2019年4月30日,沈阳农业大学生物科学技术学院植物生理与种质创新团队的阮燕晔副教授和张敖博士赴上海农业科学院作物育种栽培研究所学习考察并就玉米全基因组选择和GWAS研究在鲜食玉米种的应用进行交流研讨。


2019年4月29日-2019年4月30日,沈阳农业大学生物科学技术学院植物生理与种质创新团队的阮燕晔副教授和张敖博士赴上海农业科学院作物育种栽培研究所学习考察并就玉米全基因组选择和GWAS研究在鲜食玉米种的应用进行交流研讨。


该教程用于Visual Studio Code测试版本(无中文版选择),正式版本可以直接选择中文版。
————
VSCode是微软开发的轻量级、免费、开源、跨平台的代码编辑器,具有界面友好,功能强大等特点。不幸的是,默认情况下,该编辑器是英文的,E文不好的同学可能使用起来会有些麻烦。其实,VSCode简单3步就可以切换到中文界面。

3. 点右下角的【Restart Now】,重新打开VSCode。
4. 界面已经修改为中文。
所有样品测定完后,Stop→Ok→Save保存图像,Graph→Set axes(调整时间轴使图像在一个界面内显示),Tools→Get rate(求得图像斜率),Files→Data to ASCⅡ(保存数据)。
前面已经学过,R语言可以通过数字和运算符计算出想要的结果。但是,当我们需要对更多的数据进行计算时,比如要计算某个班级里学生的平均成绩,每次计算都输入所有班级同学的得分,然后再除以班级人数,显然看起来会很繁琐。
例子:
计算班级平均分,分数如下:
张三(95),李四(88),王二(79),这五(59),那六(92),(31),(66),(92),(73),(94),(0),(55)。
# 平均分计算
(95+88+79+59+92+31+66+92+73+94+0+55)/12
[1] 68.66667
班级人数少,计算起来还算容易,如果班级有上百人,或者人事处需要统计整个学校教师的收入情况,找出收入最高和最低的人,同时计算工资总数和平均数。那么,每次计算都输入一长串的公式显然是难以承受的。基于此,“向量”的概念被提出。该“向量”并非数学中的向量,在数学中,向量表示的是有大小有方向的量,可以用带箭头的线段来表示,线段的长度表示大小,箭头表示方向;编程中的向量一般表示一维数据的集合,它可以包含数字、字符、逻辑值甚至是三者混排,这与数学中的向量含义不同。回到上面的例子,利用向量来计算班级的平均分:
# 利用向量进行平均分计算
scores <- c(95,88,79,59,92,31,66,92,73,94,0,55) # 将一组数输入向量,用c()函数。
sum(scores)/length(scores) # sum()求和函数,length()求元素个数函数。
[1] 68.66667
mean(scores) # mean()求平均值函数。
[1] 68.66667
prod(scores) # 求乘积[1] 0
max(scores) # max()最大值函数
[1] 95
min(scores) # min()最小值函数
[1] 0
median(scores) # median()中位数函数,如果是偶数,计算中间两个数的均值。
[1] 76
quantile(scores) # quantile()求四分位数函数。
0% 25% 50% 75% 100%
0 58 76 92 95
summary(scores) # summary()总结函数,不同类型的元素结果不同。
Min. 1st Qu. Median Mean 3rd Qu. Max.
0.00 58.00 76.00 68.67 92.00 95.00
str(scores) # str()紧凑显示该对象的R语言内部结构。
num [1:12] 95 88 79 59 92 31 66 92 73 94 …
细心的同学可能会发现一个问题,向量中的考试分数并不知道是来自哪位同学的,换句话说就是不知道谁得了多少分,这对理解数据显然是不利的。R语言的设计者当然也意识到了这个问题,因此,向量是可以添加“名字”的,用到names()函数。一定要注意,给数据添加名字的顺序要和原始数据一致。
names(scores) # 直接使用names()函数,结果得到未定义NULL。
NULL
names(scores) <- c("张三","李四","王二","这五","那六") # 给向量加名字scores # 查看scores向量
张三 李四 王二 这五 那六 <NA> <NA> <NA> <NA> <NA> <NA> <NA>
95 88 79 59 92 31 66 92 73 94 0 55
(a <- c("张三"=95,"李四"=88,"王二"=79,"这五"=92)) # 直接带名字赋值向量,引号要改成英文状态,最外层加括号直接显示结果
张三 李四 王二 这五
95 88 79 92
a["张三"] # 查看张三的分数
张三
95
向量还可以进行增减和排序等操作。
<p value="<amp-fit-text layout="fixed-height" min-font-size="6" max-font-size="72" height="80">向量子集的提取同样非常有用,比如,我们要提取高于或等于九十分一个或多个界限的分数有哪几个?
length(scores) # 查看scores向量有多少个元素
[1] 12
scores[13] <- 77 # 增加第13个元素为77
scores[11] <- 20 # 将第11个元素改成20
scores <- scores[-12] # 去掉第12个元素
scores[1:3] <- c(20,20,20) # 将前3个元素改成20
scores[c(4,6,10)] # 查看第4、6、10号元素
scores[-c(4,6,10)] # 查看除了4、6、10号元素外的其他元素
scores[c(4,6,10)] <- c(40,60,80) # 将第4、6、10个元素改成40,60,80
scores[1:13] <- 100 # 1到13个元素改成100
sort(scores) # 升序排列
[1] 0 31 55 59 66 73 79 88 92 92 94 95
sort(scores,decreasing = TRUE) # 降序排列
[1] 95 94 92 92 88 79 73 66 59 55 31 0
scores[scores>=90]
[1] 95 92 92 94
scores[scores>=60]
[1] 95 88 79 92 66 92 73 94
scores[scores>=60&scores<85]
[1] 79 66 73
我相信,一定有“好事儿”的同学想验证一下,不大于60分的人和小于等于60分的子集提取方法,结果是否一样。在本例中,因为数据较少,结果一目了然,但在实际工作中,会处理大量数据,很难用肉眼分辨两组数据是否一样,这时候就用到了一个函数identical(),如果两个向量的顺序和值都相同,则得到TRUE,否则得到FALSE。
scores[scores<60]
[1] 59 31 0 55
scores[!scores>=60]
[1] 59 31 0 55
identical(scores[!scores>=60],scores[scores<60])
[1] TRUE
好了,熟练掌握上面的内容,基本上本章已经可以达标了。下面看看简单的思考题。
思考题:
9月11日,应生物科学技术学院邀请,加拿大农业与食品部渥太华研发中心(ORDC)马保罗研究员来校进行学术交流并作题为《Concurrent Improvement in Crop Productivity and Sustaining the Environment-Climate Smart Plant Nutrition》的学术报告,系统介绍精准农业提高作物生产率的整体构思,阐述了各种非生物胁迫和植物营养对作物的影响,介绍了通过反射光谱获取作物群体信息的方法、作物营养管理在养分吸收、相互作用及营养平衡等方面的进展,以及农田管理技术在提高农作物产量及环境资源高效利用中的应用等。
马保罗博士,加拿大农业和农业食品部渥太华研发中心(ORDC)作物生理学资深研究员,同时担任加拿大东部油菜规模化生产首席研究员。重点研究领域为作物生理学、栽培生理生态和营养管理。发表了170 多篇学术论文、18 篇专著章节/综述论文。马博士曾担任《Agronomy Journal》为期三年的高级副主编和《Crop Science》两届三年期的学术主编。并在2007 年被授予美国农学会学会院士,2011 年被授予加拿大农学会学会院士,2012 年被授予美国作物学会学会院士。(生物科学技术学院 张敖)

2018年8月17日-2018年8月23日,沈阳农业大学生物科学技术学院(沈阳市全基因组选择重点实验室)阮燕晔副教授和张敖博士访问了位于菲律宾马尼拉的国际水稻研究所,参加了实用单倍型作图(PHG)和卓越育种平台(EiB)的研讨会。研讨会上,康奈尔大学著名玉米统计学教授ED Bucklar及其团队详细介绍了利用单倍型作为分子标记的实现原理、优势和不足。Bucklar实验室的研究人员已经开发出将各种测序或基因分型平台获取的单核苷酸多态性标记转换为实用单倍型的工具。
卓越育种平台是国际农业磋商组织以实现发展中国家育种现代化和解决粮食及营养问题为目标而开发的开源育种信息管理分析平台,是基因组开源育种信息联盟GOBii项目的拓展。平台完成后会为公立或私营企业提供前沿的育种软件、服务和最佳事件方案。该平台包括5个模块,分别是卓越育种模块、性状发现与育种工具与服务模块、基因分型工具与服务模块、表型工具与服务模块和数据管理分析工具与服务模块。沈阳农业大学参与了由GOBii项目组领衔的基因分型、表型鉴定、质量控制、全基因组选择等方面的研究协作。


(生物科学技术学院(沈阳市全基因组选择重点实验室) 阮燕晔 张敖)
变量应该算是编程语言中应用范围最广的了,任何一个语言都离不开变量。变量,即一个容器,你将一些东西存放进去,例如:你有一个盒子,这个盒子你取名叫a,然后把6这个数字放进去,那么,不论谁拿到这个盒子,都能看见里面有一个6,当你将7放进去,同时把6取出来以后,任何人都能看到这个盒子里面有7,这就是最初的变量思想。在R语言中,变量可以存储任何东西,可以存储数值、字母、表达式、数组、数据框等等,有些概念还没有讲过,不必纠结,只要知道变量可以存储很多东西就可以了。
下面学习如何操作变量,我们设置一个变量a,然后让a变量的值为15。
a <- 15 # R语言的默认赋值时“箭头”赋值,箭头可以有方向。
15 -> a # 该结果与上面的结果一样。
a = 15 # 该结果也与上面的结果一样。
a # 查看变量内容,只需要输入变量名称即可。
当然,我们还可以将字符放入变量中。
b <- “Hello World!” # 字符串我们用英文状态下的双引号引起来,告诉程序这是字符串。
b # 查看变量b内容。
ls() # 是list的缩写,用于查看当前内存中的所有变量。
rm(a) # 当变量不再使用时,可以移除变量,用rm()函数。
rm(a,b) # 移除变量a和b,多个变量同时移除,用英文逗号隔开。
rm(list=ls()) # 移除所有变量。
在实际应用中,最好不要用单个字母来做变量名。变量名称需要易于理解,这样程序阅读或修改起来才更加容易。程序不是写一次就结束了,很多时候,程序需要根据情况不断地进行修改。
我比较习惯的命名规则为,rowNo,moneyTotal,minAlleleFrequency等。当然还有其它比较好的命名规则,如有需要请自行百度。要说明的是,不要怕变量名长,人类理解比机器理解要重要的多,随遇到满屏幕的单字母变量都会头疼,读起来会很慢。
另外,变量名称需要以字母开头,可以有下划线_或百分号%,其他符号是不允许的。
函数实际上是实现某些功能的代码块,R语言本身内置了很多函数,方便我们使用。rm是Remove的缩写,移除变量。
目前,“大刍草起源假说”在学界认可度最高。主要的原因有:
以前的史料认为,玉米起源于南美洲,直到墨西哥城考古发现八万年前的“野生玉米花粉”化石,才确立了墨西哥作为玉米发源地的地位。玉米的传播路线同样具有很大争议,学界最认可的是1492年哥伦布首次将玉米带回欧洲,再从欧洲传到全世界。
玉米进入我国主要有三种假说:
南海路传入说(先从欧洲传到菲律宾, 后由葡萄牙或菲律宾等地经商的中国人经海路传入中国)。
无论是玉米的起源还是传播路线,都缺乏有力的证据,随着各类组学、考古学和新兴科学技术的发展,希望未来能够获得一个完美的解释。
QTL IciMapping是中国农业科学院作物科学研究所开发的一款整合型的QTL作图工具。
(1)BIN,冗余标记去除;
(2)MAP,连锁图谱构建;
(3)CMP:多图谱共有标记分析;
(4)SDL,分离失真基因座图谱构建;
(5)BIP,加、显性和双基因上位效应作图;
(6)MET,执行G×E,多环境需要多个输入文件;
(7)CSL,染色体片段置换系的加、显、上位性作图;
(8)NAM,NAM群体的QTL作图。
其中,前4项功能只需要分子标记信息即可进行,后4项功能既需要分子标记信息也需要表型数据。
群体类型不同,等位基因及基因型的理论频率不同,QTL IciMapping将群体分为20个类型+一个NAM群体,如图所示:

理论频率如下表:
| Population no. | Full name | Genotype | Allele | |||
|---|---|---|---|---|---|---|
| AA | AB | BB | A | B | ||
| 1 | P1BC1F1 | 0.5 | 0.5 | 0 | 0.75 | 0.25 |
| 2 | P2BC1F1 | 0 | 0.5 | 0.5 | 0.25 | 0.75 |
| 3 | F1DH | 0.5 | 0 | 0.5 | 0.5 | 0.5 |
| 4 | F1RIL | 0.5 | 0 | 0.5 | 0.5 | 0.5 |
| 5 | P1BC1RIL | 0.75 | 0 | 0.25 | 0.75 | 0.25 |
| 6 | P2BC1RIL | 0.25 | 0 | 0.75 | 0.25 | 0.75 |
| 7 | F2 | 0.25 | 0.5 | 0.25 | 0.5 | 0.5 |
| 8 | F3 | 0.375 | 0.25 | 0.375 | 0.5 | 0.5 |
| 9 | P1BC2F1 | 0.75 | 0.25 | 0 | 0.875 | 0.125 |
| 10 | P2BC2F1 | 0 | 0.25 | 0.75 | 0.125 | 0.875 |
| 11 | P1BC2RIL | 0.875 | 0 | 0.125 | 0.875 | 0.125 |
| 12 | P2BC2RIL | 0.125 | 0 | 0.875 | 0.125 | 0.875 |
| 13 | P1BC1F2 | 0.625 | 0.25 | 0.125 | 0.75 | 0.25 |
| 14 | P2BC1F2 | 0.125 | 0.25 | 0.625 | 0.25 | 0.75 |
| 15 | P1BC2F2 | 0.8125 | 0.125 | 0.0625 | 0.875 | 0.125 |
| 16 | P2BC2F2 | 0.0625 | 0.125 | 0.8125 | 0.125 | 0.875 |
| 17 | P1BC1DH | 0.75 | 0 | 0.25 | 0.75 | 0.25 |
| 18 | P2BC1DH | 0.25 | 0 | 0.75 | 0.25 | 0.75 |
| 19 | P1BC2DH | 0.875 | 0 | 0.125 | 0.875 | 0.125 |
| 20 | P2BC2DH | 0.125 | 0 | 0.875 | 0.125 | 0.875 |
AA和BB是两个亲本贡献的基因型。
| Genotype | AA | AB | BB | Missing | AB + AA | AB + BB |
|---|---|---|---|---|---|---|
| Encoding options | 2 | 1 | 0 | − 1 | 12, D | 10, R |
| A | H | B | X | AH, HA | BH, HB | |
| AA | AB | BB | XX | AX, XA | BX, XB | |
| BA | * | A*, *A | B*, *B | |||
| ** | A_, _A | B_, _B | ||||
| Encoding by number | 2 | 1 | 0 | − 1 | 12 | 10 |
| Encoding by single letter | A | H | B | X | D | R |
| Encoding by double letters | AA | AB | BB | XX | AX | XB |
比较推荐使用数字编码规则,即亲本1=2,亲本2=0,杂合=1,缺失=-1。在孟德尔遗传条件下,无法区分杂合和纯和显性的情况下,使用12和10,SNP标记由于可以确定和某一个亲本一致,因此用不到这两项。
| MTP | 在遗传群体中,不同标记和分离失真的卡方检验的频率 |
| SUM | 构建组合连锁图谱的概要 |
| LOD | 组合连锁关系(针对独立继承)的成对标记LOD值测试的下三角矩阵 |
| REC | 标记成对组合重组率估计估计的下三角矩阵 |
| DIS | 从组合重组率计算标记成对遗传距离(cM)的下三角矩阵 |
| SDL | SDL功能的输入文件,用于分离失真基因座作图 |
| BIP | BIL功能的输入文件,可在添加表型文件后进行加性和上位性作图 |
| MET | MET功能的输入文件,添加多环境数据后,计算QTL与环境互作 |
染色体片段置换系(CSS,chromosome segment substitution),也叫做渐渗系(introgression lines),通常需要经过反复回交(有时使用分子标记来辅助选择),突出目标片段,而使非有效基因座上的基因替换为已知背景。在理想状态下,每个片段中包含一个QTL,去除背景,很容易验证是否片段中携带QTL。然而,CSS群体的构建需要很长时间,而且工作量巨大。在低世代的CSS中,有些片段含有多个QTL。在高强度的选择下,基因和标记的频率不遵循孟德尔遗传定律。因此,CSS的QTL作图方法采用基于逐步回归的似然比检测(stepwise regression-based likelihood ratio test,RSTEP-LRT)。
NAM(nested association mapping)群体,是巢式关联作图群体,由一个共同亲本和多个亲本产生的若干RIL家系组成的多重交叉组配设计群体。对NAM群体的QTL方法叫联合复合区间作图群体( joint inclusive composite interval mapping,JICIM)。