TWASandGWAS中GBS filtering and GWAS（1）

F:\文章代碼\TWASandGWAS\GBS filtering and GWAS

README.TXT

請檢查幻燈片“Vitamaize_update_Gorelab_Ames_GBS_filtering_20191122.pptx”中關于阿姆斯（Ames）ID處理流程的詳細信息。

文件夾“Ames_ID_processing”包含了用于處理阿姆斯ID的文件和R腳本。
文件“Ames_GBS_raw_filtering_for_imputation.txt”包含了用于過濾GBS（基因組簡化測序）基因型數據以供填補使用的Linux命令。

Ames_ID_processing解讀

1. 導入和準備數據

在這部分，腳本首先加載了處理數據所需的所有R庫。這些庫提供了數據操作、數據可視化、數據重塑等功能。接著，腳本讀取了包含Pedigree信息、位置信息等的CSV文件。這些文件包含了進行后續分析所需的基礎數據。

目的：準備分析環境，導入原始數據。

2. Pedigree信息比較

這一部分的目的是比較不同來源的Pedigree信息，確保數據的一致性。腳本讀取了由Maria和Laura提供的Pedigree信息，并進行了比較，找出兩者之間的差異。

目的：確保Pedigree信息的一致性，為后續分析打下基礎。

3. 處理位置信息

在這一部分，腳本處理了一個包含位置信息的文件。這個文件包含了樣本的位置和來源信息，這些信息對于后續的分析非常重要。

目的：整理和清洗位置信息，確保每個樣本的位置信息準確無誤。

4. 連接Maria端和Laura端的行

這一部分的目的是根據位置信息，將Maria端的樣本與Laura端的源ID進行匹配。這樣可以幫助研究者理解樣本的來源和背景。

目的：將樣本與其來源ID進行關聯，為后續的遺傳分析提供信息。

5. 過濾和處理GBS記錄

在這部分，腳本從Panzea_ZeaGBSv2.7中檢索完整的GBS記錄列表，并根據需要移除甜玉米和爆裂玉米的記錄。

目的：篩選出感興趣的樣本，排除不需要的樣本，如甜玉米和爆裂玉米。

6. 處理和生成鎖定列表

這一部分涉及移除特定的行（如"ae"行和七個通道），生成一個鎖定列表，這個列表包含了經過篩選的樣本。

目的：生成一個清潔的樣本列表，用于后續的遺傳分析。

7. 檢索共同行的GBS記錄

在這一部分，腳本從Panzea_ZeaGBSv2.7中檢索共同行的GBS記錄列表。

目的：獲取用于后續分析的共同樣本的GBS記錄。

8. 數據清洗和驗證

這一部分的目的是驗證和清洗數據，確保數據的質量。這包括修正NPGS中的不一致性，確保數據的一致性。

目的：確保數據的準確性和一致性，為后續分析提供高質量的數據。

9. 生成最終的GBS記錄列表

最后，腳本生成最終的GBS記錄列表，這個列表包含了所有經過篩選和驗證的樣本，可以用于后續的遺傳分析。

目的：生成一個最終的樣本列表，用于后續的遺傳分析。

每個部分都是數據處理流程中不可或缺的一環，確保了從原始數據到分析就緒數據的轉換過程的準確性和可靠性。如果你對某個部分有更具體的問題或需要更詳細的解釋，請告訴我，我可以進一步解釋。

1. 導入和準備數據

首先，我們需要加載R中處理數據所需的庫。這些庫提供了數據操作、數據可視化、數據重塑等功能。

# 加載必要的庫
library(gdata)
library(reshape2)
library(ggplot2)
library(Hmisc)
library(dplyr)
library(gridExtra)

解釋：

gdata：用于數據處理。
reshape2：用于數據重塑，例如從寬格式轉換為長格式。
ggplot2：一個流行的數據可視化庫，用于創建高質量的圖形。
Hmisc：提供額外的數據操作工具。
dplyr：提供數據操作函數，特別是對數據框（data frame）的操作。
gridExtra：提供額外的圖形布局功能。

這些庫是進行數據分析和可視化的基礎，確保了后續步驟可以順利進行。

接下來，我們將讀取包含Pedigree信息、位置信息等的CSV文件。這些文件包含了進行后續分析所需的基礎數據。

# 讀取Pedigree信息文件
coder <- read.csv('entry_coder_isu.csv',header = TRUE)
coder_ex <- coder[which(coder$Pedigree!= "B73"),]
coder_p_dup <- coder_ex[duplicated2(coder_ex$Pedigree),] # 12 pedigree_dup
coder_s_dup <- coder_ex[duplicated2(coder_ex$Source),] # zero source_dup
write.csv(coder_p_dup ,file="coder_p_dup.csv", row.names = F)

解釋：

read.csv：讀取CSV文件，header = TRUE表示文件的第一行包含列名。
coder <- coder[which(coder$Pedigree!= "B73"),]：篩選出Pedigree不是"B73"的行。
duplicated2：用于找出重復值的行。
coder_p_dup <- coder_ex[duplicated2(coder_ex$Pedigree),]：找出Pedigree列中重復的行。
write.csv：將結果寫入新的CSV文件。

這部分代碼的目的是從Pedigree信息中篩選出非"B73"的樣本，并找出其中的重復項，然后將這些重復項寫入新的文件中。

2. Pedigree信息比較

這一部分的目的是確保不同來源的Pedigree信息的一致性，并處理不同年份的數據。

### d2015_vb
d2015_vb <- read.csv('Ames 2015 - BVits Final_021518.csv',header = TRUE)  # 1802 lines
d2015_vb <- d2015_vb[,c(6,8,9)]
d2015_vb$Range_Pass_2015_vb <- paste(d2015_vb$Range,d2015_vb$Pass,sep="_")
d2015_vb_c <- merge(d2015_vb,coder,by.x="Range_Pass_2015_vb",by.y="Range_Pass_2015")
d2015_vb_c_dif <- d2015_vb_c[which(as.character(d2015_vb_c$Pedigree.x)!=as.character(d2015_vb_c$Pedigree.y)),]### d2015_ve
d2015_ve <- read.csv('AMES15 - TOCOCHROMATICS - REVISED_022219.csv',header = TRUE) # 1801 lines
d2015_ve <- d2015_ve[,c(6,8,9)]
d2015_ve$Range_Pass_2015_ve <- paste(d2015_ve$Range,d2015_ve$Pass,sep="_")
d2015_ve_c <- merge(d2015_ve,coder,by.x="Range_Pass_2015_ve",by.y="Range_Pass_2015")
d2015_ve_c_dif <- d2015_ve_c[which(as.character(d2015_ve_c$Pedigree.x)!=as.character(d2015_ve_c$Pedigree.y)),]d2015_dif <-merge(d2015_vb_c_dif,d2015_ve_c_dif,by.x="Range_Pass_2015_vb",by.y="Range_Pass_2015_ve",all=T)
write.csv(d2015_dif ,file="d2015_dif.csv", row.names = F)

解釋：

read.csv：讀取CSV文件，header = TRUE表示文件的第一行包含列名。
d2015_vb <- d2015_vb[,c(6,8,9)]：選擇特定的列（第6、8、9列）。
paste：將Range和Pass列合并為一個新列Range_Pass_2015_vb。
merge：合并兩個數據框（data frame），基于共同的列。
which(as.character(d2015_vb_c$Pedigree.x)!=as.character(d2015_vb_c$Pedigree.y))：找出兩個Pedigree列值不匹配的行。
write.csv：將結果寫入新的CSV文件。

這部分代碼首先讀取2015年的兩個數據集（d2015_vb和d2015_ve），然后分別與coder數據框合并，以比較Pedigree信息。接著，找出兩個Pedigree列值不匹配的行，并將這些行寫入新的CSV文件d2015_dif.csv中。

3. 處理位置信息

這一部分的目的是檢查和處理包含樣本位置信息的文件，以確保后續分析中樣本位置的準確性。

### 讀取原始文件
key <- read.csv('location_to_Source_15-17_key.csv', header = TRUE) # 3790 rows
### 添加 "id" 列
key$id <- paste(key$Year, key$Range, key$Pass, sep="_")
key_dup <- key[duplicated(key$id),] # 10 rows with duplications
write.csv(key_dup, file="key_dup.csv", row.names = F)
### 移除10行重復項
key_u <- distinct(key, id, .keep_all = TRUE) # 3790 - 10 = 3780 rows
write.csv(key_u, file="location_to_Source_15-17_key_update.csv", row.names = F)### 讀取更新后的key文件
key <- read.csv('location_to_Source_15-17_key_update.csv', header = TRUE) # 3780 rows

解釋：

read.csv：讀取包含位置信息的CSV文件。
paste：創建一個新的列id，該列是通過將Year、Range和Pass列合并生成的，用于唯一標識每一行。
duplicated：找出id列中重復的行。
write.csv：將包含重復項的數據框寫入新的CSV文件key_dup.csv中。
distinct：移除數據框中的重復行，保留唯一的行。
key_u：更新后的無重復項數據框。

這部分代碼首先讀取位置信息文件，然后創建一個新的唯一標識列id。接著，代碼找出并移除了重復的行，最后將更新后的數據框寫入新的CSV文件中。

4. 連接Maria端的行與Laura端的源ID

這一部分的目的是將Maria端的行數據與Laura端的源ID進行匹配，以便將樣本的基因型信息與其來源聯系起來。

### 加載必要的庫
library(gdata)
library(reshape2)
library(ggplot2)
library(Hmisc)
library(dplyr)
library(gridExtra)### 讀取更新后的key文件
key <- read.csv('location_to_Source_15-17_key_update.csv', header = TRUE) # 3780 rows### d2015_vb
d2015_vb <- read.csv('Ames 2015 - BVits Final_021518.csv', header = TRUE)  # 1802 lines
d2015_vb <- d2015_vb[,c(6,8,9)]
d2015_vb$id_2015_vb <- paste("2015", d2015_vb$Range, d2015_vb$Pass, sep="_")
d2015_vb_k <- merge(d2015_vb, key, by.x="id_2015_vb", by.y="id")
d2015_vb_k_dif <- d2015_vb_k[which(as.character(d2015_vb_k$Pedigree.x) != as.character(d2015_vb_k$Pedigree.y)),]### d2015_ve
d2015_ve <- read.csv('AMES15 - TOCOCHROMATICS - REVISED_022219.csv', header = TRUE) # 1801 lines
d2015_ve <- d2015_ve[,c(6,8,9)]
d2015_ve$id_2015_ve <- paste("2015", d2015_ve$Range, d2015_ve$Pass, sep="_")
d2015_ve_k <- merge(d2015_ve, key, by.x="id_2015_ve", by.y="id")
d2015_ve_k_dif <- d2015_ve_k[which(as.character(d2015_ve_k$Pedigree.x) != as.character(d2015_ve_k$Pedigree.y)),]### d2017_vb
d2017_vb <- read.csv('Ames 2017_BVits_FINAL_d3Thiamine corrected and plate info corrected_022219.csv', header = TRUE) # 1748 lines
d2017_vb <- d2017_vb[,c(6,8,9)]
d2017_vb$id_2017_vb <- paste("2017", d2017_vb$Range, d2017_vb$Pass, sep="_")
d2017_vb_k <- merge(d2017_vb, key, by.x="id_2017_vb", by.y="id")
d2017_vb_k_dif <- d2017_vb_k[which(as.character(d2017_vb_k$Pedigree.x) != as.character(d2017_vb_k$Pedigree.y)),]### d2017_ve
d2017_ve <- read.csv('AMES 2017_TOCOCHROMATICS_FINAL_.csv', header = TRUE) # 1738 lines
d2017_ve <- d2017_ve[,c(6,8,9)]
d2017_ve$id_2017_ve <- paste("2017", d2017_ve$Range, d2017_ve$Pass, sep="_")
d2017_ve_k <- merge(d2017_ve, key, by.x="id_2017_ve", by.y="id")
d2017_ve_k_dif <- d2017_ve_k[which(as.character(d2017_ve_k$Pedigree.x) != as.character(d2017_ve_k$Pedigree.y)),]### 比較Pedigree信息，完成！

解釋：

首先，代碼加載了必要的R庫，以便進行數據處理和分析。
然后，代碼讀取了更新后的key文件，該文件包含了樣本的位置和來源信息。
對于2015年和2017年的數據，代碼分別讀取了兩個數據集（d2015_vb和d2017_vb，d2015_ve和d2017_ve），并選擇了特定的列。
接著，代碼為每個數據集創建了一個新的id列，該列是通過合并Range和Pass列生成的，用于唯一標識每一行。
然后，代碼將每個數據集與key文件合并，以將樣本的基因型信息與其來源信息關聯起來。
最后，代碼找出了兩個Pedigree列值不匹配的行，并進行了處理。

5. 生成Ames ID完整列表以用于GBS提取

這一部分的目的是生成一個完整的Ames ID列表，用于后續的基因組簡化測序（GBS）提取工作。

### 加載必要的庫
library(gdata)
library(reshape2)
library(ggplot2)
library(Hmisc)
library(dplyr)
library(gridExtra)### 讀取更新后的key文件
key <- read.csv('location_to_Source_2015-17_key_update.csv', header = TRUE) # 3780 rows### d2015_vb
d2015_vb <- read.csv('Ames_2015_vb_adding_ID.csv', header = TRUE)  # 1802 lines
d2015_vb$id_2015_vb <- paste("2015", d2015_vb$Range, d2015_vb$Pass, sep="_")
d2015_vb_k <- merge(d2015_vb, key, by.x="id_2015_vb", by.y="id")
d2015_vb_k2 <- d2015_vb_k[,c(2:6,31:33,27:28,8:26)]
colnames(d2015_vb_k2)[5] <- "Source"
colnames(d2015_vb_k2)[7] <- "Pedigree"
colnames(d2015_vb_k2)[8] <- "ID"
colnames(d2015_vb_k2)[12] <- "Range"
colnames(d2015_vb_k2)[13] <- "Pass"
### 檢查每個行是否只被測量了一次
d2015_vb_k2_dup <- d2015_vb_k2[duplicated2(d2015_vb_k2$PedId),] # yes
write.csv(d2015_vb_k2 ,file="Ames_2015_vb_adding_ID.csv", row.names = F)#### d2015_ve
d2015_ve <- read.csv('Ames_2015_ve_adding_ID.csv', header = TRUE) # 1801 lines
d2015_ve$id_2015_ve <- paste("2015", d2015_ve$Range, d2015_ve$Pass, sep="_")
d2015_ve_k <- merge(d2015_ve, key, by.x="id_2015_ve", by.y="id")
d2015_ve_k2 <- d2015_ve_k[,c(2:6,31:33,27:28,8:26)]
colnames(d2015_ve_k2)[5] <- "Source"
colnames(d2015_ve_k2)[7] <- "Pedigree"
colnames(d2015_ve_k2)[8] <- "ID"
colnames(d2015_ve_k2)[12] <- "Range"
colnames(d2015_ve_k2)[13] <- "Pass"
### 檢查每個行是否只被測量了一次
d2015_ve_k2_dup <- d2015_ve_k2[duplicated2(d2015_ve_k2$PedId),] # yes
write.csv(d2015_ve_k2 ,file="Ames_2015_ve_adding_ID.csv", row.names = F)#### d2017_vb
d2017_vb <- read.csv('Ames_2017_vb_adding_ID.csv', header = TRUE) # 1748 lines
d2017_vb$id_2017_vb <- paste("2017", d2017_vb$Range, d2017_vb$Pass, sep="_")
d2017_vb_k <- merge(d2017_vb, key, by.x="id_2017_vb", by.y="id")
d2017_vb_k2 <- d2017_vb_k[,c(2:6,31:33,27:28,8:26)]
colnames(d2017_vb_k2)[5] <- "Source"
colnames(d2017_vb_k2)[7] <- "Pedigree"
colnames(d2017_vb_k2)[8] <- "ID"
colnames(d2017_vb_k2)[12] <- "Range"
colnames(d2017_vb_k2)[13] <- "Pass"
### 檢查每個行是否只被測量了一次
d2017_vb_k2_dup <- d2017_vb_k2[duplicated2(d2017_vb_k2$PedId),] # 10 samples were measured two times
### 為了保持一致性，移除那些在Plate 19上測量兩次的10個樣本
remove_id <- d2017_vb_k2_dup[which(d2017_vb_k2_dup$LCMS.Plate=="Ames17-BV-19L"),c(1:3)]
all_row <- NULL # 獲取行號
for (i in 1:10){tmp_row <- which(d2017_vb_k2$LCMS.Plate=="Ames17-BV-19L" & d2017_vb_k2$PedId==remove_id$PedId[i] )all_row <- c(all_row,tmp_row)
}
d2017_vb_k3 <- d2017_vb_k2[-all_row,]
d2017_vb_k3_dup <- d2017_vb_k3[duplicated2(d2017_vb_k3$PedId),]#### d2017_ve
d2017_ve <- read.csv('Ames_2017_ve_adding_ID.csv', header = TRUE) # 1738 lines
d2017_ve$id_2017_ve <- paste("2017", d2017_ve$Range, d2017_ve$Pass, sep="_")
d2017_ve_k <- merge(d2017_ve, key, by.x="id_2017_ve", by.y="id")
d2017_ve_k2 <- d2017_ve_k[,c(2:6,30:32,26:27,8:25)]
colnames(d2017_ve_k2)[5] <- "Source"
colnames(d2017_ve_k2)[7] <- "Pedigree"
colnames(d2017_ve_k2)[8] <- "ID"
colnames(d2017_ve_k2)[12] <- "Range"
colnames(d2017_ve_k2)[13] <- "Pass"
### 檢查每個行是否只被測量了一次
d2017_ve_k2_dup <- d2017_ve_k2[duplicated2(d2017_ve_k2$PedId),] # all samples were measured one time### 連接Maria端的行與Laura端的源ID，完成！

解釋：

首先，代碼加載了必要的R庫，以便進行數據處理和分析。
然后，代碼讀取更新后的key文件，該文件包含了樣本的位置和來源信息。
對于2015年和2017年的數據，代碼分別讀取了兩個數據集（d2015_vb和d2017_vb，d2015_ve和d2017_ve），并選擇了特定的列。
接著，代碼為每個數據集創建了一個新的id列，該列是通過合并Range和Pass列生成的，用于唯一標識每一行。
然后，代碼將每個數據集與key文件合并，以將樣本的基因型信息與其來源信息關聯起來。
最后，代碼檢查了每個行是否只被測量了一次，并進行了相應的處理。

這部分代碼的目的是將樣本的基因型信息與其來源信息進行匹配，以便進行后續的遺傳分析。如果你對這部分代碼理解清楚了，我們可以繼續講解下一部分。

6. 從Panzea_ZeaGBSv2.7中檢索完整的GBS記錄列表

這一部分的目的是從Panzea_ZeaGBSv2.7數據集中檢索與Ames ID相關的所有GBS記錄，以便進行后續的基因型分析。

### 加載必要的庫
library(gdata)
library(reshape2)
library(ggplot2)
library(Hmisc)
library(dplyr)
library(gridExtra)### 讀取Ames ID完整列表
ames_isu <- read.csv(file="ames_isu_id.csv",stringsAsFactors=FALSE)### 讀取Panzea_ZeaGBSv2.7的記錄列表
v2.7 <- read.csv(file="Panzea_ZeaGBSv2.7_id.csv")### 清洗和格式化ID列
ames_isu$ID_clean <- ames_isu$ID
ames_isu$ID_clean <- gsub(' ','',ames_isu$ID_clean)
ames_isu$ID_clean <- tolower(ames_isu$ID_clean)v2.7$acc_v2.7_clean <- v2.7$acc_v2.7
v2.7$acc_v2.7_clean <- gsub(' ','',v2.7$acc_v2.7_clean)
v2.7$acc_v2.7_clean <- gsub('-','',v2.7$acc_v2.7_clean)
v2.7$acc_v2.7_clean <- gsub('_','',v2.7$acc_v2.7_clean)
v2.7$acc_v2.7_clean <- gsub('\\(','',v2.7$acc_v2.7_clean)
v2.7$acc_v2.7_clean <- gsub('\\)','',v2.7$acc_v2.7_clean)
v2.7$acc_v2.7_clean <- tolower(v2.7$acc_v2.7_clean)### 合并數據框
ames_isu_v2.7 <- merge(ames_isu, v2.7, by.x="ID_clean", by.y="acc_v2.7_clean",all = T)### 過濾掉沒有GBS記錄的樣本
ames_isu_v2.7 <- ames_isu_v2.7[which(!is.na(ames_isu_v2.7$INDV)),]### 保存結果
write.csv(ames_isu_v2.7 ,file="ames_isu_v2.7_all.csv", row.names = F)

解釋：

首先，代碼加載了必要的R庫。
然后，代碼讀取了Ames ID的完整列表（ames_isu）和Panzea_ZeaGBSv2.7的記錄列表（v2.7）。
接著，代碼對ID列進行了清洗和格式化，以確保它們可以正確匹配。
之后，代碼合并了這兩個數據框，以將Ames ID與對應的GBS記錄關聯起來。
最后，代碼過濾掉了沒有GBS記錄的樣本，并將結果保存到CSV文件中。

這部分代碼的目的是生成一個包含所有Ames ID及其對應GBS記錄的完整列表，這對于后續的遺傳分析非常重要。如果你對這部分代碼理解清楚了，我們可以繼續講解下一部分。

7. 移除甜玉米和爆裂玉米

這一部分的目的是從GBS記錄列表中移除甜玉米和爆裂玉米的記錄，因為這些樣本可能對分析結果產生干擾。

### 加載必要的庫
library(gdata)
library(reshape2)
library(ggplot2)
library(Hmisc)
library(dplyr)
library(gridExtra)### 讀取Ames ID完整列表
ames_isu <- read.csv(file="ames_isu_id.csv", header = TRUE)### 讀取包含玉米品種信息的文件
ames_2013gb <- read.csv('13059_2013_103_MOESM1_ESM.csv', header = TRUE) # 品種列表### 為Ames ID添加清洗后的格式
ames_isu$ID_cleanformat <- ames_isu$ID
ames_isu$ID_cleanformat <- gsub(' ','', ames_isu$ID_cleanformat)
ames_isu$ID_cleanformat <- tolower(ames_isu$ID_cleanformat)### 合并品種列表和Ames ID列表
ames_isu_gb <- merge(ames_isu, ames_2013gb, by.x="ID_cleanformat", by.y="Accesion.N", all=T)### 修正Cinta列表中的重復錯誤
ames_isu_gb_dup <- ames_isu_gb[duplicated(ames_isu_gb$ID_cleanformat),]
ames_isu_gb <- ames_isu_gb[-which(ames_isu_gb$ID_cleanformat == "ames27101" & ames_isu_gb$Breeding.program == "Ontario"),]
ames_isu_gb <- ames_isu_gb[-which(ames_isu_gb$ID_cleanformat == "pi543850" & ames_isu_gb$Breeding.program == "Other"),]
ames_isu_gb <- ames_isu_gb[,-c(5:11)]
ames_isu_gb$Comments_merged <- paste(ames_isu_gb$Comments, ames_isu_gb$Pop.structure, sep="_")
table(ames_isu_gb$Comments_merged)
ames_isu_gb$Comments_isu_gb <- "other"
ames_isu_gb$Comments_note <- "ok"
ames_isu_gb$Comments_isu_gb[which(ames_isu_gb$Comments_merged == "popcorn or sweet corn_unclassified")] <- "popcorn or sweet corn_unclassified"
ames_isu_gb$Comments_note[which(ames_isu_gb$Comments_merged == "popcorn or sweet corn_unclassified")] <- "in_question"
ames_isu_gb$Comments_isu_gb[which(ames_isu_gb$Comments_merged == "popcorn_popcorn")] <- "popcorn"
ames_isu_gb$Comments_note[which(ames_isu_gb$Comments_merged == "popcorn_popcorn")] <- "in_question"
ames_isu_gb$Comments_isu_gb[which(ames_isu_gb$Comments_merged == "popcorn_unclassified")] <- "popcorn_unclassified"
ames_isu_gb$Comments_note[which(ames_isu_gb$Comments_merged == "popcorn_unclassified")] <- "in_question"
ames_isu_gb$Comments_isu_gb[which(ames_isu_gb$Comments_merged == "sweet corn_NA")] <- "sweet corn_NA"
ames_isu_gb$Comments_note[which(ames_isu_gb$Comments_merged == "sweet corn_NA")] <- "in_question"
ames_isu_gb$Comments_isu_gb[which(ames_isu_gb$Comments_merged == "sweet corn_non-stiff stalk")] <- "sweet corn_non-stiff stalk"
ames_isu_gb$Comments_note[which(ames_isu_gb$Comments_merged == "sweet corn_non-stiff stalk")] <- "in_question"
ames_isu_gb$Comments_isu_gb[which(ames_isu_gb$Comments_merged == "sweet corn_stiff stalk")] <- "sweet corn_stiff stalk"
ames_isu_gb$Comments_note[which(ames_isu_gb$Comments_merged == "sweet corn_stiff stalk")] <- "in_question"
ames_isu_gb$Comments_isu_gb[which(ames_isu_gb$Comments_merged == "sweet corn_sweet corn")] <- "sweet corn"
ames_isu_gb$Comments_isu_gb[which(ames_isu_gb$Comments_merged == "sweet corn_unclassified")] <- "sweet corn_unclassified"
ames_isu_gb$Comments_note[which(ames_isu_gb$Comments_merged == "sweet corn_unclassified")] <- "in_question"ames_isu_gb <- ames_isu_gb[,-c(4:6)]
in_ques_38 <- ames_isu_gb[which(ames_isu_gb$Comments_note=="in_question"),] ### 38 accessions
write.csv(in_ques_38,file="ames_isu_in_question_38.csv", row.names = F)### 檢查這些訪問號在NPGS中的信息并修正差異
npgs_38 <- read.csv('ames_isu_in_question_38_npgs.csv', header = TRUE,stringsAsFactors=FALSE) ### 36 of the 38 in_question were fixed
ames_isu_gb_ok <- ames_isu_gb[which(ames_isu_gb$Comments_note=="ok"),]for (i in c(1:38)) {if (npgs_38$NPGS[i]=="sweet corn"|npgs_38$NPGS[i]=="popcorn") { npgs_38$Comments_note[i] <- "ok" npgs_38$Comments_isu_gb[i] <- npgs_38$NPGS[i]}if (npgs_38$Comments_note[i]=="in_question") { npgs_38$Comments_note[i] <- paste(npgs_38$Comments_note[i],"_NPGS_", npgs_38$NPGS[i],sep="")}
}ames_isu_gb_npgs38 <- rbind(ames_isu_gb_ok,npgs_38[,c(1:5)])
write.csv(ames_isu_gb_npgs38,file="ames_isu_gb_npgs38.csv", row.names = F)### 合并Di的列表信息
ames_isu_gb_npgs38 <- read.csv('ames_isu_gb_npgs38.csv', header = TRUE)
check_83 <- read.csv('lines_to_check_sweet_pop_vitamaize_master.csv', header = TRUE)
ames_isu_gb_npgs38 <- merge(ames_isu_gb_npgs38, check_83, by="ID_cleanformat",all=T)
ames_isu_gb_npgs38 <- ames_isu_gb_npgs38[which(!is.na(ames_isu_gb_npgs38$ID)),]
table(ames_isu_gb_npgs38$GRIN)for (i in c(1:1762)) {if (!is.na(ames_isu_gb_npgs38$GRIN[i]) & ames_isu_gb_npgs38$GRIN[i]=="popcorn") { ames_isu_gb_npgs38$Comments_isu_gb[i] <- ames_isu_gb_npgs38$GRIN[i]ames_isu_gb_npgs38$Comments_note[i] <- "ok" }if (!is.na(ames_isu_gb_npgs38$GRIN[i]) & ames_isu_gb_npgs38$GRIN[i]=="unknown") { ames_isu_gb_npgs38$Comments_note[i] <- "in_question_NPGS_unknown"}
}
write.csv(ames_isu_gb_npgs38[,c(1:5)],file="ames_isu_gb_npgs38_npgs8.csv", row.names = F) ### with four accessions still in question

解釋：

代碼首先讀取了Ames ID完整列表和包含玉米品種信息的文件。
接著，代碼對ID列進行了清洗和格式化，以確保它們可以正確匹配。
然后，代碼合并了這兩個數據框，以將Ames ID與對應的品種信息關聯起來。
代碼修正了Cinta列表中的重復錯誤，移除了重復的行。
代碼檢查了這些訪問號在NPGS中的信息并修正了差異。
最后，代碼合并了Di的列表信息，生成了最終的品種列表。

這部分代碼的目的是從GBS記錄列表中移除甜玉米和爆裂玉米的記錄，以確保分析結果的準確性。如果你對這部分代碼理解清楚了，我們可以繼續講解下一部分。

8. 比較甜玉米和爆裂玉米列表與Laura的列表

這一部分的目的是確保我們移除的甜玉米和爆裂玉米列表與Laura提供的列表一致，從而保證數據的一致性。

### 加載必要的庫
library(gdata)
library(reshape2)
library(ggplot2)
library(Hmisc)
library(dplyr)
library(gridExtra)### 讀取Ames ID完整列表
ames_isu <- read.csv('ames_isu_gb_npgs38_npgs83.csv', header = TRUE)### 讀取Laura提供的甜玉米和爆裂玉米列表
check_83 <- read.csv('sweet.and.pop.for.Di.Source_Laura.csv', header = TRUE)### 從我們的列表中移除Laura列表中的樣本
ames_isu_gb <- ames_isu_gb[!(ames_isu_gb$ID %in% check_83$Source),]### 比較我們的列表和Laura的列表，確保一致性
ames_isu_gb$ID_cleanformat <- gsub(' ','', ames_isu_gb$ID)
ames_isu_gb$ID_cleanformat <- tolower(ames_isu_gb$ID_cleanformat)check_83$ID_cleanformat <- gsub(' ','', check_83$Source)
check_83$ID_cleanformat <- tolower(check_83$ID_cleanformat)### 合并兩個列表，找出差異
mer <- merge(ames_isu_gb, check_83, by="ID_cleanformat", all=T)### 檢查是否有不一致的樣本
mer_oth <- mer[is.na(mer$Gen),]
mer_oth_sp <- mer_oth[which(mer_oth$Comments_isu_gb!="other"),] ### sp did not contain additional sweet or pop### 檢查一致性
mer_con <- mer[!is.na(mer$Gen),]
mer_con$comment_merge <- paste(mer_con$Comments_isu_gb,"#",mer_con$Comments,sep="")
table(mer_con$comment_merge)
ques <- mer_con[which(mer_con$Comments_isu_gb=="other"),] ### these 2 lines are in-question from the Comments_note### 完成：我們的甜玉米和爆裂玉米列表與Laura的相同

解釋：

首先，代碼讀取了我們之前生成的Ames ID列表（ames_isu_gb_npgs38_npgs3.csv）和Laura提供的甜玉米和爆裂玉米列表（sweet.and.pop.for.Di.Source_Laura.csv）。
接著，代碼從我們的列表中移除了Laura列表中的樣本，確保我們的列表與Laura的列表一致。
然后，代碼對ID列進行了清洗和格式化，以確保它們可以正確匹配。
之后，代碼合并了兩個列表，以找出任何可能的差異。
最后，代碼檢查了合并后的列表，找出不一致的樣本，并確認我們的列表與Laura的列表是否一致。

這部分代碼的目的是確保我們的甜玉米和爆裂玉米列表與Laura提供的列表一致，從而保證數據的一致性和可靠性。如果你對這部分代碼理解清楚了，我們可以繼續講解下一部分。