Introduction

In this vignette, we will implement the genomic selection in a Maize pipeline. Two scenarios will be implemented after a fifteen-year Burn-In period. The first scenario will account only for truncated phenotypic selection. In the second scenario, we will implement genomic selection in the pipeline. So, in years 3 and 4, the individuals will be selected based on the estimated breeding values that came from a genomic model (ridge regression model). In addition, the parents will be selected based on ‘ebv’ and came from year 3` inbreds.

Package and files for recording the outputs

# Packages
require(AlphaSimR)
## Loading required package: AlphaSimR
## Loading required package: R6
# Creating the files to record the results
MeanG_pop = matrix(NA, 35)
MeanA_pop = matrix(NA, 35)
VarG_pop  = matrix(NA, 35)
Accuracy  = matrix(NA, 35)

Create parents and fill the pipeline

# Parameters files
source("aux_files/GlobalParameters.R")

# Create the parents
source("aux_files/Create_parents.R")

# FillPipeline
source("aux_files/FillPipeline.R")
## FillPipeline year: 1 of 6
## FillPipeline year: 2 of 6
## FillPipeline year: 3 of 6
## FillPipeline year: 4 of 6
## FillPipeline year: 5 of 6
## FillPipeline year: 6 of 6
# Environmental covariate
P = runif(burninYears+futureYears)

Burn-In period

for(year in 1:burninYears){
p = P[year]
cat("Working on year:",year,"\n")
source("aux_files/UpdateParents.R")  # Pick new parents based on last year's data
source("aux_files/UpdateTesters.R")  # Pick new testers and hybrid parents
source("aux_files/AdvanceYearPS.R")  # Advances yield trials by a year
source("aux_files/WriteRecordsGS.R") # Write records for GS predictions
source("aux_files/UpdateResults.R")  # Track summary data

}
## Working on year: 1 
##  Writing records
## Loading required package: AGHmatrix
## Working on year: 2 
##  Writing records 
## Working on year: 3 
##  Writing records 
## Working on year: 4 
##  Writing records 
## Working on year: 5 
##  Writing records 
## Working on year: 6 
##  Writing records 
## Working on year: 7 
##  Writing records 
## Working on year: 8 
##  Writing records 
## Working on year: 9 
##  Writing records 
## Working on year: 10 
##  Writing records 
## Working on year: 11 
##  Writing records 
## Working on year: 12 
##  Writing records 
## Working on year: 13 
##  Writing records 
## Working on year: 14 
##  Writing records 
## Working on year: 15 
##  Writing records
# Saving the info

save.image("results/BURNIN.RData")

Phenotypic selection program

# 3.0 Loading the scenarios
load("results/BURNIN.RData")


for(year in (burninYears+1):(burninYears+futureYears)){
p = P[year]

# 3.1 Loop
cat("Working on year:",year,"\n")
source("aux_files/UpdateParents.R")  # Pick new parents based on last year's data
source("aux_files/UpdateTesters.R")  # Pick new testers and hybrid parents
source("aux_files/AdvanceYearPS.R")  # Advances yield trials by a year
source("aux_files/UpdateResults.R")  # Track summary data
  
}
## Working on year: 16 
## Working on year: 17 
## Working on year: 18 
## Working on year: 19 
## Working on year: 20 
## Working on year: 21 
## Working on year: 22 
## Working on year: 23 
## Working on year: 24 
## Working on year: 25 
## Working on year: 26 
## Working on year: 27 
## Working on year: 28 
## Working on year: 29 
## Working on year: 30 
## Working on year: 31 
## Working on year: 32 
## Working on year: 33 
## Working on year: 34 
## Working on year: 35
# 3.2 Recording results
output1 = data.frame(scenario=rep("PS", 35), # Scenario name
                     MeanG_pop,
                     Accuracy,
                     MeanA_pop,
                     VarG_pop,
                     stringsAsFactors=FALSE)

# 3.3 Saving the results as RDS
saveRDS(output1,"results/Results_PS.rds")

Genomic selection program

# 3.0 Loading the scenarios
load("results/BURNIN.RData")


for(year in (burninYears+1):(burninYears+futureYears)){
p = P[year]

cat("Working on year:",year,"\n")
if(year == (burninYears+1)){
  source('aux_files/fillGS.R')
}
source("aux_files/UpdateParents_GS.R")  # Pick new parents based on last year's data
source("aux_files/UpdateTesters.R")  # Pick new testers and hybrid parents
source("aux_files/AdvanceYearGS.R")  # advance the populations
source("aux_files/WriteRecordsGS.R")  # Track records for GS
source("aux_files/UpdateResults.R")  # Track summary data
  
}
## Working on year: 16 
##  Writing records 
## Working on year: 17 
##  Writing records 
## Working on year: 18 
##  Writing records 
## Working on year: 19 
##  Writing records 
## Working on year: 20 
##  Writing records 
## Working on year: 21 
##  Writing records 
## Working on year: 22 
##  Writing records 
## Working on year: 23 
##  Writing records 
## Working on year: 24 
##  Writing records 
## Working on year: 25 
##  Writing records 
## Working on year: 26 
##  Writing records 
## Working on year: 27 
##  Writing records 
## Working on year: 28 
##  Writing records 
## Working on year: 29 
##  Writing records 
## Working on year: 30 
##  Writing records 
## Working on year: 31 
##  Writing records 
## Working on year: 32 
##  Writing records 
## Working on year: 33 
##  Writing records 
## Working on year: 34 
##  Writing records 
## Working on year: 35 
##  Writing records
# 3.2 Recording results
output1 = data.frame(scenario=rep("GS", 35), # Scenario name
                     MeanG_pop,
                     Accuracy,
                     MeanA_pop,
                     VarG_pop,
                     stringsAsFactors=FALSE)

# 3.3 Saving the results as RDS
saveRDS(output1,"results/Results_GS.rds")

Plotting the results

# Loading the results
Scenario1 = readRDS("results/Results_GS.rds")
Scenario2 = readRDS("results/Results_PS.rds")


# Plot hybrid mean of genetic values over time
meanRanges = range(c(Scenario1$MeanG_pop[15:35], Scenario2$MeanG_pop[15:35]))
plot(x = 15:35, y = Scenario1$MeanG_pop[15:35], type = "l", col = "black", lwd = 3,
     xlab = "Year", ylab = "Mean of genetic values", ylim = meanRanges)
lines(x = 15:35, y = Scenario2$MeanG_pop[15:35], type = "l", col = "blue", lwd = 3)
legend(x = "topleft", legend = c('GS', 'Pheno'), title = "Scenarios",
       lwd = 3, lty = c(1, 1), col = c("black", "blue"), bty = "n")

References

  1. University of Florida, ↩︎

  2. University of Florida, ↩︎