10 min read

常规遗传评估需要的数据格式和整理步骤

文件格式

  • 为了确保在Win和Linux系统下均可使用,确保数据和系谱保存格式为utf-8。
  • 育种数据分析,主要用到2个文件,系谱文件和数据文件。如果开展基因组选择,会使用基于分子标记信息的G或者H矩阵或者逆矩阵文件。在ASReml中,通常以稀疏矩阵格式提供(三列)逆矩阵文件。
  • 系谱文件和数据文件,通常保存为.csv文本格式。csv格式就是逗号分隔的文本格式。具体保存方法,在excel中输入数据后,点击保存,选择csv格式。通常在中文Win下excel生成的csv文件是GBK格式,需要用Ultraedit或者Notepad++重新保存为utf-8格式。
  • 数据文件名表示方式:性状前缀+世代名称Data.csv,譬如中国对虾2015年生长测试数据文件,可以用g2015G11Data.csv表示。
  • 系谱文件名表示方式:性状前缀++世代名称Ped.csv,譬如中国对虾2015年系谱文件,可以用g2015G11Ped.csv表示。
  • 如果是基于分子亲缘关系逆矩阵分析,通常扩展名为giv文件,每一列通过空格或者制表符分割。

通用规则

  • 数据按照行保存,每一行是一个个体的记录。每一列表示不同的变量,如个体编号,父本编号,性别,出生年份等信息。
  • 第一行一般是表头,表示该列字段的意义。譬如第一列是个体编号,通常在第一行用AnimalID来表示。第二列是父本编号,通常用SireID表示。
  • 约定俗成的规定,表头一般用英文字符如AnimalID等来表示,因为某些统计分析或者数据处理软件软件不兼容中文。
  • 数据文件中,缺失的数据元素,用NA表示。譬如,如果看不出某只虾的性别,用NA代替。如果某只虾的性别鉴定不出来,也用NA表示。
  • 数据文件中,性状关联字段,如收获体重标记时体重等数值字段,务必请检查数据范围,修正明显超出范围的错误值。譬如收获体重均值为10g,如果某只虾的体重为300g、1000g,明显是输入错误;如果存在0.01g,很大概率也是输入错误。
  • 对于性状、测试池等分类字段,可以利用Excel中的筛选功能,查看每个字段有多少个水平。譬如性别列,如果出现Male、male、 Female、F、female等,需要统一修正为仅包括两种性别。
  • 当用Excel整理数据时,一定要确保单元格中不能为空,如果没有值,用NA填充。系谱文件例外,见下文。
  • 如果用ASReml软件进行遗传评估,生成的扩展名为.ass的数据汇总文件务必要首先查看,核对数据是否准确。
  • 推荐利用R包DataExplorer包对待分析数据进行检查。

系谱文件

  • 系谱文件有三个字段,分别为个体编号、父本编号和母本编号;
  • 系谱文件个体数据要按照规则排序:一个个体的父母本系谱信息必须出现在该个体前边;
  • 如果个体的父本和母本未知,用0代替。父本和母本要么都为0,如果只知道已授精的母本,父本也要编一个号码;
  • 如果系谱文件和数据文件共用一个文件,那么文件的前三列将是个体编号、父本编号和母本编号;
  • 如果父本编号和母本编号两列中的个体没有在个体编号列中出现,将会默认为奠基者个体。
  • 在第四列加入个体的性别是一个好的建议。
  • 可以利用optiSel、visPedigree等R包实现对系谱文件的整理
  • visPedigree包的中文安装和使用说明

文件保存

在数据分析过程中会产生大量的文件,为了保证数据分析的可靠性,请务必使用以下文件夹层次保存文件:

PWS2010G01 #育种项目顶层文件夹
—–Data
——–Input #保存原始的系谱和数据文件
——–Output #保存R脚本清洗后的数据和系谱等文件
—–ASReml #保存as文件和ASReml生成的输出文件
—–R #保存清洗数据的R脚本文件
—–BLUPF90 #生成H矩阵的BLUPF90脚本

如果使用ASReml分析数据,as文件中务必使用!OUTFOLDER限定符将不同Path的内容输出到单独文件夹中。关于R文件的编码格式,请参考这里

数据处理实例

knitr::opts_chunk$set(echo = TRUE)
rm(list=ls()) #清空已存在的各种对象和变量
library(data.table) 
library(visPedigree)
library(DataExplorer)

1.读取系谱和数据文件

数据文件系谱文件可以点击下载。
运行R代码前的注意事项和Rstudio编辑器配置和运行等,参见另外一个blog
下载后,跟R代码放在同一目录下,如果不打算修改下述代码,那么把数据和系谱文件放在一个新建目录datasets下。datasets文件夹和R脚本文件在同一个目录下。

gPed_dt <-
  fread(
    file = "datasets/g2016G0to2017G1Ped.csv",
    sep = ",",
    stringsAsFactors = FALSE,
    header = TRUE
  )
gData_dt <-
  fread(
    file = "datasets/g2016G0to2017G1Data.csv",
    sep = ",",
    stringsAsFactors = FALSE,
    header = TRUE
  )
gPed_dt$Sex[gPed_dt$Sex == "0"] <- NA

2. 数据检查和整理

2.1 数据检查

检查数据中分类以及连续性变量的分布,主要是发现是否有异常值。

根据柱状图,可以查看世代、池塘、性别等分类变量的水平数。如果性别有三个,很容易可以查看出来。

根据直方图,可以检查体重、体长、日龄等连续变量的分布,如果有异常值,也非常容易查看。

连续性变量之间的散点图和相关图,用来发现不同变量间的趋势,也是非常有用的工具。

#检查分类变量
plot_bar(gData_dt)
## 4 columns ignored with more than 50 categories.
## AnimalID: 11413 categories
## SireID: 90 categories
## DamID: 120 categories
## FamilyID: 120 categories

#检查连续性变量
plot_histogram(gData_dt)

#检查连续性变量与收获体重间的散点图
plot_scatterplot(split_columns(gData_dt)$continuous, by = "M2BW")

#检查连续性变量间的相关性
plot_correlation(gData_dt, type = "continuous")

2.2 数据清洗

数据分析的需求是经常变动和碎片化的,因此需要掌握熟练的数据清洗技能。数据清洗主流的R包有dplyrdata.table等包。网络上关于二者的比较也有不少。前者的优势是有一个较为完善的数据处理体系,具体参见tidyverse包合集。data.table的优势是特别适合处理大数据,读入、处理和写入的速度非常快。考虑到测序和分型数据都是上G的数据量,推荐利用data.table包。data.table包理解起来,一开始可能有些门槛,不过跟着教程,很快就会熟悉了。我也跟着官方英文教程写了几个学习的blog,供大家参考。百度、知乎上也有不少data.table的使用方法介绍,大家自行参考,具体data.table的使用就不在此赘述了。这里举几个在育种数据分析中用到data.table的例子:

2.2.1 数据提取

(1)提取雄虾数据
有时候,审稿人可能会要求我们单独分析雄虾数据的遗传参数,这就要求我们提取雄虾的数据,利用data.table,实现起来很简单。注意%chin%操作符是data.table单独针对字符类型重构的,速度较快。数值类型,可以直接用>、<等符号。

gData_dt[SexID %chin% c("Female")]
##              AnimalID          SireID           DamID Generation
##    1: OXHMB17G010D158 OXHMB16G000A389 OXHMB16G000B754        G01
##    2: OXHMB17G010D261 OXHMB16G000A012 MRTZD16G000F031        G01
##    3: OXHMB17G0101515 OXHMB16G000F297 OXHMB16G000E364        G01
##    4: OXHMB17G0109465 OXHMB16G000A343 MRTZD16G000F026        G01
##    5: OXHMB17G0103104 OXHMB16G000A765  OXHMB16G00X110        G01
##   ---                                                           
## 4477: OXHMB17G0108813 OXHMB16G000A275 MRTZD16G000F001        G01
## 4478: OXHMB17G0000146 MRTZD16B000M002 MRTZD16B000F002        G01
## 4479: OXHMB17G0110114 MRTZD16G000M052 MRTZD16G000F052        G01
## 4480: OXHMB17G0106486 OXHMB16G000F132 OXHMB16G000F371        G01
## 4481: OXHMB17G010F177 MRTZD16G000M021 OXHMB16G000B013        G01
##              FamilyID TankID  SexID  Birthday    S2Date    M1Date M1BW
##    1: OXHMB17G0100482    P11 Female 4/22/2117 5/14/2117  7/9/2117 16.4
##    2: OXHMB17G0100891    P11 Female 4/24/2117 5/15/2117 7/10/2117 14.7
##    3: OXHMB17G0100211    P10 Female 4/19/2117  5/9/2117  7/8/2117 15.2
##    4: OXHMB17G0100881    P10 Female 4/18/2117 5/11/2117  7/8/2117 13.6
##    5: OXHMB17G0100323    P10 Female  5/1/2117 5/23/2117 7/12/2117 13.0
##   ---                                                                 
## 4477: OXHMB17G0100831    P10 Female 4/20/2117  5/9/2117  7/8/2117 16.4
## 4478: OXHMB17G0001021    P10 Female 5/10/2117 5/29/2117 7/14/2117  8.7
## 4479: OXHMB17G0100933    P10 Female 4/30/2117 5/22/2117 7/12/2117 11.3
## 4480: OXHMB17G0100611    P10 Female 4/28/2117 5/20/2117 7/11/2117 12.8
## 4481: OXHMB17G0100771    P10 Female 4/14/2117  5/4/2117  7/6/2117 15.4
##          M2Date M2Age1 M2Age2 M2BW M2BL
##    1: 9/28/2117    159    137  430  113
##    2: 9/28/2117    157    136  488  117
##    3: 10/9/2117    173    153  337  111
##    4: 10/9/2117    174    151  273  103
##    5: 10/9/2117    161    139  355  109
##   ---                                  
## 4477: 10/9/2117    172    153  362  111
## 4478: 10/9/2117    152    133  276  102
## 4479: 10/9/2117    162    140  157   85
## 4480: 10/9/2117    164    142  201   98
## 4481: 10/9/2117    178    158  689  133

举一反三,如何提取体重大于700的个体记录集?

gData_dt[M2BW >700]
##           AnimalID          SireID           DamID Generation        FamilyID
## 1: OXHMB17G010D260 OXHMB16G000F302 OXHMB16G000B783        G01 OXHMB17G0100461
## 2: OXHMB17G010D342 OXHMB16G000A316 OXHMB16G000F439        G01 OXHMB17G0100131
## 3: OXHMB17G010D872 MRTZD16G000M021 OXHMB16G000B013        G01 OXHMB17G0100771
## 4: OXHMB17G010F039 OXHMB16G000C612 OXHMB16G000B679        G01 OXHMB17G0100232
## 5: OXHMB17G010F063 MRTZD16G000M001  OXHMB16G00X112        G01 OXHMB17G0100731
## 6: OXHMB17G010F087 MRTZD16G000M021 OXHMB16G000B013        G01 OXHMB17G0100771
## 7: OXHMB17G010F216 OXHMB16G000A805 OXHMB16G000F485        G01 OXHMB17G0100441
## 8: OXHMB17G010F254 OXHMB16G000C612 OXHMB16G000B679        G01 OXHMB17G0100232
## 9: OXHMB17G010F899 OXHMB16G000C612 OXHMB16G000B679        G01 OXHMB17G0100232
##    TankID SexID  Birthday    S2Date    M1Date M1BW    M2Date M2Age1 M2Age2 M2BW
## 1:    P11  Male 4/14/2117  5/4/2117  7/7/2117 19.7 9/28/2117    167    147  774
## 2:    P11  Male 4/13/2117  5/1/2117  7/6/2117 17.8 9/28/2117    168    150  706
## 3:    P10  Male 4/14/2117  5/4/2117  7/6/2117 15.4 10/9/2117    178    158  725
## 4:    P10  Male 4/14/2117  5/2/2117  7/6/2117 17.4 10/9/2117    178    160  744
## 5:    P10  Male 4/25/2117 5/16/2117 7/11/2117 13.8 10/9/2117    167    146  741
## 6:    P10  Male 4/14/2117  5/4/2117  7/6/2117 15.4 10/9/2117    178    158  703
## 7:    P10  Male 4/15/2117  5/6/2117  7/7/2117 16.1 10/9/2117    177    156  770
## 8:    P10  Male 4/14/2117  5/2/2117  7/6/2117 17.4 10/9/2117    178    160  850
## 9:    P10  Male 4/14/2117  5/2/2117  7/6/2117 17.4 10/9/2117    178    160  709
##    M2BL
## 1:  130
## 2:  125
## 3:  130
## 4:  130
## 5:  130
## 6:  128
## 7:  126
## 8:  132
## 9:  129

(2)提取指定家系的雄虾数据
如果想提取特定几个家系,如编号为OXHMB17G0100482、OXHMB17G0100891、OXHMB17G0100211家系的雄虾,如何处理?&是R中的逻辑运算符,表示两边条件都满足时,才返回TRUE。

gData_dt[(FamilyID %chin% c("OXHMB17G0100482", "OXHMB17G0100891", "OXHMB17G0100211")) & (SexID %chin% c("Male"))]
##             AnimalID          SireID           DamID Generation        FamilyID
##   1: OXHMB17G0105286 OXHMB16G000A389 OXHMB16G000B754        G01 OXHMB17G0100482
##   2: OXHMB17G0101510 OXHMB16G000F297 OXHMB16G000E364        G01 OXHMB17G0100211
##   3: OXHMB17G010A413 OXHMB16G000A012 MRTZD16G000F031        G01 OXHMB17G0100891
##   4: OXHMB17G010A418 OXHMB16G000F297 OXHMB16G000E364        G01 OXHMB17G0100211
##   5: OXHMB17G010A438 OXHMB16G000A012 MRTZD16G000F031        G01 OXHMB17G0100891
##  ---                                                                           
## 191: OXHMB17G0109666 OXHMB16G000A012 MRTZD16G000F031        G01 OXHMB17G0100891
## 192: OXHMB17G0109667 OXHMB16G000A012 MRTZD16G000F031        G01 OXHMB17G0100891
## 193: OXHMB17G0109668 OXHMB16G000A012 MRTZD16G000F031        G01 OXHMB17G0100891
## 194: OXHMB17G0109669 OXHMB16G000A012 MRTZD16G000F031        G01 OXHMB17G0100891
## 195: OXHMB17G0109670 OXHMB16G000A012 MRTZD16G000F031        G01 OXHMB17G0100891
##      TankID SexID  Birthday    S2Date    M1Date M1BW    M2Date M2Age1 M2Age2
##   1:    P11  Male 4/22/2117 5/14/2117  7/9/2117 16.4 9/28/2117    159    137
##   2:    P11  Male 4/19/2117  5/9/2117  7/8/2117 16.2 9/28/2117    162    142
##   3:    P11  Male 4/24/2117 5/15/2117 7/10/2117 14.7 9/28/2117    157    136
##   4:    P11  Male 4/19/2117  5/9/2117  7/8/2117 16.2 9/28/2117    162    142
##   5:    P11  Male 4/24/2117 5/15/2117 7/10/2117 14.7 9/28/2117    157    136
##  ---                                                                        
## 191:    P10  Male 4/24/2117 5/15/2117 7/10/2117 15.4 10/9/2117    168    147
## 192:    P10  Male 4/24/2117 5/15/2117 7/10/2117 15.4 10/9/2117    168    147
## 193:    P10  Male 4/24/2117 5/15/2117 7/10/2117 15.4 10/9/2117    168    147
## 194:    P10  Male 4/24/2117 5/15/2117 7/10/2117 15.4 10/9/2117    168    147
## 195:    P10  Male 4/24/2117 5/15/2117 7/10/2117 15.4 10/9/2117    168    147
##      M2BW M2BL
##   1:  169   89
##   2:  417  111
##   3:  522  118
##   4:  604  121
##   5:  515  120
##  ---          
## 191:  333  103
## 192:  360  109
## 193:  241  102
## 194:  414  109
## 195:  444  113

(3)插入新列
譬如给整个数据集加入一个流水号,用于后续数据整理。请注意,:=会直接在原数据集中插入新的一列。如果你不想改变原有数据集,需要提前拷贝原数据集,再插入新列。

gData_new_dt <- copy(gData_dt)
gData_new_dt[,":="(SerialNum=.I)]
gData_new_dt
##               AnimalID          SireID           DamID Generation
##     1: OXHMB17G010A338 OXHMB16G000F302 OXHMB16G000B783        G01
##     2: OXHMB17G010A339 OXHMB16G000C364  OXHMB16G00X751        G01
##     3: OXHMB17G010A340 OXHMB16G000C863 OXHMB16G000B274        G01
##     4: OXHMB17G0101354 OXHMB16G000C625 OXHMB16G000F625        G01
##     5: OXHMB17G010A341 OXHMB16G000A300 MRTZD16G000F021        G01
##    ---                                                           
## 11409: OXHMB17G0104360 OXHMB16G000C514 OXHMB16G000E734        G01
## 11410: OXHMB17G0110468 MRTZD16G000M071 MRTZD16G000F071        G01
## 11411: OXHMB17G0000260 MRTZD16B000M003 MRTZD16B000F003        G01
## 11412: OXHMB17G0000513 MRTZD16B000M005 MRTZD16B000F005        G01
## 11413: OXHMB17G0104470 OXHMB16G000F703 OXHMB16G000F694        G01
##               FamilyID TankID SexID  Birthday    S2Date    M1Date M1BW
##     1: OXHMB17G0100461    P11  Male 4/14/2117  5/4/2117  7/7/2117 19.7
##     2: OXHMB17G0100411    P11  Male 4/14/2117  5/1/2117  7/6/2117 16.0
##     3: OXHMB17G0100661    P11  Male 4/24/2117 5/15/2117 7/10/2117 15.2
##     4: OXHMB17G0100201    P11  Male 4/21/2117 5/14/2117  7/9/2117 15.6
##     5: OXHMB17G0100871    P11  Male 4/12/2117  5/1/2117  7/6/2117 18.2
##    ---                                                                
## 11409: OXHMB17G0100422    P10  Male  5/1/2117 5/23/2117 7/12/2117 14.2
## 11410: OXHMB17G0100971    P10  Male 4/25/2117 5/16/2117 7/11/2117 15.1
## 11411: OXHMB17G0001031    P10  Male 5/10/2117 5/29/2117 7/14/2117  9.1
## 11412: OXHMB17G0001051    P10  Male 5/10/2117 5/30/2117 7/14/2117  7.7
## 11413: OXHMB17G0100433    P10  Male  5/8/2117 5/28/2117 7/13/2117 12.0
##           M2Date M2Age1 M2Age2 M2BW M2BL SerialNum
##     1: 9/28/2117    167    147  395  113         1
##     2: 9/28/2117    167    150  408  112         2
##     3: 9/28/2117    157    136  465  120         3
##     4: 9/28/2117    160    137  355  112         4
##     5: 9/28/2117    169    150  581  123         5
##    ---                                            
## 11409: 10/9/2117    161    139  214   94     11409
## 11410: 10/9/2117    167    146  156   85     11410
## 11411: 10/9/2117    152    133  142   83     11411
## 11412: 10/9/2117    152    132  174   89     11412
## 11413: 10/9/2117    154    134   50   59     11413

.I是data.table中定义的特殊字符,有特殊的含义,这里表示会生成一个从1开始,间隔为1的序列。因此也可以这样实现:

gData_new_dt[,":="(SerialNum=seq(from=1, to=nrow(gData_new_dt), by=1))]
gData_new_dt
##               AnimalID          SireID           DamID Generation
##     1: OXHMB17G010A338 OXHMB16G000F302 OXHMB16G000B783        G01
##     2: OXHMB17G010A339 OXHMB16G000C364  OXHMB16G00X751        G01
##     3: OXHMB17G010A340 OXHMB16G000C863 OXHMB16G000B274        G01
##     4: OXHMB17G0101354 OXHMB16G000C625 OXHMB16G000F625        G01
##     5: OXHMB17G010A341 OXHMB16G000A300 MRTZD16G000F021        G01
##    ---                                                           
## 11409: OXHMB17G0104360 OXHMB16G000C514 OXHMB16G000E734        G01
## 11410: OXHMB17G0110468 MRTZD16G000M071 MRTZD16G000F071        G01
## 11411: OXHMB17G0000260 MRTZD16B000M003 MRTZD16B000F003        G01
## 11412: OXHMB17G0000513 MRTZD16B000M005 MRTZD16B000F005        G01
## 11413: OXHMB17G0104470 OXHMB16G000F703 OXHMB16G000F694        G01
##               FamilyID TankID SexID  Birthday    S2Date    M1Date M1BW
##     1: OXHMB17G0100461    P11  Male 4/14/2117  5/4/2117  7/7/2117 19.7
##     2: OXHMB17G0100411    P11  Male 4/14/2117  5/1/2117  7/6/2117 16.0
##     3: OXHMB17G0100661    P11  Male 4/24/2117 5/15/2117 7/10/2117 15.2
##     4: OXHMB17G0100201    P11  Male 4/21/2117 5/14/2117  7/9/2117 15.6
##     5: OXHMB17G0100871    P11  Male 4/12/2117  5/1/2117  7/6/2117 18.2
##    ---                                                                
## 11409: OXHMB17G0100422    P10  Male  5/1/2117 5/23/2117 7/12/2117 14.2
## 11410: OXHMB17G0100971    P10  Male 4/25/2117 5/16/2117 7/11/2117 15.1
## 11411: OXHMB17G0001031    P10  Male 5/10/2117 5/29/2117 7/14/2117  9.1
## 11412: OXHMB17G0001051    P10  Male 5/10/2117 5/30/2117 7/14/2117  7.7
## 11413: OXHMB17G0100433    P10  Male  5/8/2117 5/28/2117 7/13/2117 12.0
##           M2Date M2Age1 M2Age2 M2BW M2BL SerialNum
##     1: 9/28/2117    167    147  395  113         1
##     2: 9/28/2117    167    150  408  112         2
##     3: 9/28/2117    157    136  465  120         3
##     4: 9/28/2117    160    137  355  112         4
##     5: 9/28/2117    169    150  581  123         5
##    ---                                            
## 11409: 10/9/2117    161    139  214   94     11409
## 11410: 10/9/2117    167    146  156   85     11410
## 11411: 10/9/2117    152    133  142   83     11411
## 11412: 10/9/2117    152    132  174   89     11412
## 11413: 10/9/2117    154    134   50   59     11413

2.2.2 数据统计

(1)获得雌雄虾体重的均值

gData_dt[, .(MeanM2BW=mean(M2BW)), by=.(SexID)]
##     SexID MeanM2BW
## 1:   Male 365.4177
## 2: Female 295.4953

(2)进一步统计雌雄虾在不同池塘间的体重均值:

gData_dt[, .(MeanM2BW=mean(M2BW)), by=.(SexID, TankID)]
##     SexID TankID MeanM2BW
## 1:   Male    P11 350.4550
## 2: Female    P11 459.0000
## 3: Female    P10 295.4223
## 4:   Male    P10 376.7736

(3)继续增加需求,在统计变量中,加入体长:

gData_dt[, .(MeanM2BW=mean(M2BW), MeanM2BL=mean(M2BL)), by=.(SexID, TankID)]
##     SexID TankID MeanM2BW MeanM2BL
## 1:   Male    P11 350.4550 105.5232
## 2: Female    P11 459.0000 115.0000
## 3: Female    P10 295.4223 103.3449
## 4:   Male    P10 376.7736 108.2390

(3)需求一直在变:进一步统计P11池塘收获了多少尾母虾?

gData_dt[, .(MeanM2BW=mean(M2BW), MeanM2BL=mean(M2BL), N=.N), by=.(TankID, SexID)]
##    TankID  SexID MeanM2BW MeanM2BL    N
## 1:    P11   Male 350.4550 105.5232 2991
## 2:    P11 Female 459.0000 115.0000    2
## 3:    P10 Female 295.4223 103.3449 4479
## 4:    P10   Male 376.7736 108.2390 3941

.N也是data.table中定义的一个特殊符号,表示每个分组的行数。

(4)思考如何满足这样一个需求?提取的数据集,要求每个家系每个池塘母虾数量大于30尾。

gData_dt[, ":="(Size=.N), by=.(FamilyID, TankID, SexID)]
setorder(gData_dt,-Size)
gData_dt[Size >= 30]
##              AnimalID          SireID           DamID Generation
##    1: OXHMB17G0100507 OXHMB16G000A414  OXHMB16G00X015        G01
##    2: OXHMB17G010E349 OXHMB16G000A414  OXHMB16G00X015        G01
##    3: OXHMB17G0100508 OXHMB16G000A414  OXHMB16G00X015        G01
##    4: OXHMB17G0100509 OXHMB16G000A414  OXHMB16G00X015        G01
##    5: OXHMB17G010E394 OXHMB16G000A414  OXHMB16G00X015        G01
##   ---                                                           
## 7903: OXHMB17G0105768 OXHMB16G000F063 OXHMB16G000E752        G01
## 7904: OXHMB17G0102731 OXHMB16G000A116 OXHMB16G000E255        G01
## 7905: OXHMB17G0106033 OXHMB16G000C751 OXHMB16G000E570        G01
## 7906: OXHMB17G0102732 OXHMB16G000A116 OXHMB16G000E255        G01
## 7907: OXHMB17G0102733 OXHMB16G000A116 OXHMB16G000E255        G01
##              FamilyID TankID  SexID  Birthday    S2Date    M1Date M1BW
##    1: OXHMB17G0100072    P10 Female 4/23/2117 5/15/2117 7/10/2117 17.3
##    2: OXHMB17G0100072    P10 Female 4/23/2117 5/15/2117 7/10/2117 17.3
##    3: OXHMB17G0100072    P10 Female 4/23/2117 5/15/2117 7/10/2117 17.3
##    4: OXHMB17G0100072    P10 Female 4/23/2117 5/15/2117 7/10/2117 17.3
##    5: OXHMB17G0100072    P10 Female 4/23/2117 5/15/2117 7/10/2117 17.3
##   ---                                                                 
## 7903: OXHMB17G0100522    P10   Male 4/21/2117 5/11/2117  7/9/2117 13.9
## 7904: OXHMB17G0100302    P10   Male 4/24/2117 5/18/2117 7/11/2117 13.7
## 7905: OXHMB17G0100561    P10   Male 4/23/2117 5/15/2117 7/10/2117 11.1
## 7906: OXHMB17G0100302    P10   Male 4/24/2117 5/18/2117 7/11/2117 13.7
## 7907: OXHMB17G0100302    P10   Male 4/24/2117 5/18/2117 7/11/2117 13.7
##          M2Date M2Age1 M2Age2 M2BW M2BL Size
##    1: 10/9/2117    169    147  294  104   76
##    2: 10/9/2117    169    147  328  107   76
##    3: 10/9/2117    169    147  284  103   76
##    4: 10/9/2117    169    147  306  105   76
##    5: 10/9/2117    169    147  354  109   76
##   ---                                       
## 7903: 10/9/2117    171    151  509  118   30
## 7904: 10/9/2117    168    144  285   98   30
## 7905: 10/9/2117    169    147  395  112   30
## 7906: 10/9/2117    168    144  508  120   30
## 7907: 10/9/2117    168    144  420  115   30

3.预处理系谱、提纯系谱、数字化系谱

visPedigree包中的tidyped()函数,主要有两个功能:

  • 包括了自动补全不在系谱ID列中的个体号,以及提纯功能。
  • 把与当前数据集中个体没有亲缘关系的系谱个体,全部剔除掉
#预处理系谱,主要是把父母本的系谱设置为0
gPed_tidy_dt <- tidyped(ped = gPed_dt, cand  = gData_dt$AnimalID)
## Warning in checkped(ped, addgen): In the sire or dam column, Blank, Zero,
## asterisk, or character NA is recognized as a missing parent and is replaced with
## the missing value NA.
gPed_tidy_dt
##                    Ind            Sire             Dam    Sex Offspring Year
##     1: OLHBB14B000F001            <NA>            <NA> female      TRUE 2014
##     2: OLHBB14B000F002            <NA>            <NA> female      TRUE 2014
##     3: OLHBB14B000F003            <NA>            <NA> female      TRUE 2014
##     4: OLHBB14B000F004            <NA>            <NA> female      TRUE 2014
##     5: OLHBB14B000F005            <NA>            <NA> female      TRUE 2014
##    ---                                                                      
## 11914: OXHMB17G0110820 OLHBC16G000M088 OLHBC16G000F088   male     FALSE 2017
## 11915: OXHMB17G0110821 OLHBC16G000M088 OLHBC16G000F088   male     FALSE 2017
## 11916: OXHMB17G0110822 OLHBC16G000M088 OLHBC16G000F088   male     FALSE 2017
## 11917: OXHMB17G0110823 OLHBC16G000M088 OLHBC16G000F088   male     FALSE 2017
## 11918: OXHMB17G0110824 OLHBC16G000M088 OLHBC16G000F088   male     FALSE 2017
##          Breed  Cand Gen IndNum SireNum DamNum
##     1:      DP FALSE   1      1       0      0
##     2:      DP FALSE   1      2       0      0
##     3:      DP FALSE   1      3       0      0
##     4:      DP FALSE   1      4       0      0
##     5:      DP FALSE   1      5       0      0
##    ---                                        
## 11914: unknown  TRUE   4  11914     350    332
## 11915: unknown  TRUE   4  11915     350    332
## 11916: unknown  TRUE   4  11916     350    332
## 11917: unknown  TRUE   4  11917     350    332
## 11918: unknown  TRUE   4  11918     350    332

整理后系谱的包含的个体数,要少于原始系谱。

4.画出个体的系谱树

侧面验证一下提纯后的系谱,是否正确。

从数据集中最后一尾虾,画出系谱图。

这里需要注意的是,visped()函数要求个体的父本和母本,要么双亲独有,要么都没有,不能一个有,一个没有。

实际上可能会存在这样一个情况,像两只已经交尾的母虾,有可能来自同一尾公虾,但是即便是这种情况,也需要给父本一个随机编号,假定他们是各不相同的。

即便是设置为缺失,应该也是这种情况。

one_ped_dt <-
  tidyped(ped = gPed_tidy_dt, cand = gData_dt$AnimalID[nrow(gData_dt)])
visped(ped = one_ped_dt)

5.输出系谱

输出系谱时,设置缺失的父本和母本编号为0。

fwrite()函数来自data.table包。

gPed_tidy_dt[is.na(Sire), Sire := "0"]
gPed_tidy_dt[is.na(Dam), Dam := "0"]
fwrite(
  gPed_tidy_dt,
  file = "datasets/g2016G0to2017G1Pedpruned.csv",
  sep = ",",
  row.names = F,
  na = "NA",
  quote = FALSE
)