Reading teosinte data

library(tidyverse)
library(cowplot)
library(qtl)

teosinte.xo<-read.table(file="xo_ZeaGBSv27raw_RareAllelesC2TeoCurated_depth_AGPv4_filtered0210.txt",header=TRUE,stringsAsFactors=FALSE)
head(teosinte.xo)
##   taxon chr  start     end
## 1 A0151   1 337863  369777
## 2 A0415   1 337863  369777
## 3 A1190   1 337863  369777
## 4 A2718   1 337863  369777
## 5 A0134   1 337863 1035713
## 6 A0557   1 334927 1035713
teosinte.parentage<-read.table(file="latest_parentage_info/C2_Teosinte_Parentage_0210.txt",header=TRUE,stringsAsFactors=FALSE,sep="\t")
head(teosinte.parentage)
##    Plot        Mom        Dad    Cross
## 1 A0003 PC_N14_ID2 PC_N14_ID2     self
## 2 A0004 PC_J48_ID2 PC_J48_ID2     self
## 3 A0008 PC_J13_ID1 PC_K60_ID1 outcross
## 4 A0013 PC_N10_ID2 PC_N10_ID2     self
## 5 A0014 PC_K60_ID1 PC_K60_ID1     self
## 6 A0015 PC_J10_ID1 PC_M15_ID1 outcross
names(teosinte.parentage)<-c("Taxa","Mother","Father","cross")

Transforming teosinte data

teosinte.merged.xo<-merge(teosinte.xo,teosinte.parentage,by.x="taxon",by.y="Taxa")


teosinte.child.counts<- teosinte.parentage %>% 
  group_by(Mother) %>%
  summarise(family_progeny_count = length(Mother))

teosinte.selfed.counts<-teosinte.parentage[teosinte.parentage$Mother==teosinte.parentage$Father,] %>%
   group_by(Mother) %>%
  summarise(family_selfed_count = length(Mother))

teosinte.xo.by.mom<- teosinte.merged.xo %>% 
  group_by(Mother) %>%
  summarise(family_xo_count = length(Mother))
 

 
teosinte.merged.xo<-inner_join(teosinte.merged.xo,teosinte.child.counts)

teosinte.merged.xo<-inner_join(teosinte.merged.xo,teosinte.selfed.counts)

teosinte.merged.xo<-inner_join(teosinte.merged.xo,teosinte.xo.by.mom)

teosinte.merged.moms<-inner_join(teosinte.xo.by.mom,teosinte.child.counts)

teosinte.merged.moms<-inner_join(teosinte.merged.moms,teosinte.selfed.counts)
head(teosinte.merged.moms)
## # A tibble: 6 × 4
##       Mother family_xo_count family_progeny_count family_selfed_count
##        <chr>           <int>                <int>               <int>
## 1 PC_I05_ID1            2791                  105                  23
## 2 PC_I06_ID1            1793                   75                  38
## 3 PC_I08_ID1            2241                   84                  40
## 4 PC_I11_ID2            3381                  129                  44
## 5 PC_I50_ID1            2431                   92                  13
## 6 PC_I50_ID2            2366                   73                  55
dim(teosinte.merged.moms)
## [1] 52  4
teosinte.xo.over.progeny.count<-data.frame(rep("teosinte",dim(teosinte.merged.moms)[1]),teosinte.merged.moms$family_xo_count/teosinte.merged.moms$family_progeny_count)
names(teosinte.xo.over.progeny.count)<-c("teo.or.maize","xo_over_family_count")

head(teosinte.xo.over.progeny.count)
##   teo.or.maize xo_over_family_count
## 1     teosinte             26.58095
## 2     teosinte             23.90667
## 3     teosinte             26.67857
## 4     teosinte             26.20930
## 5     teosinte             26.42391
## 6     teosinte             32.41096

Reading maize landraces data

landraces.parentage<-read.table(file="latest_parentage_info/2017_feb19_landrace_parentage.txt",header=TRUE,stringsAsFactors=FALSE,sep="\t")
head(landraces.parentage)
##               plant_ID                                      Mother
## 1 13FL990001:250304331 171_11_Tuxpeno_El_Aguacatito_Mex_:250358103
## 2 13FL990002:250304332  167_5_Tuxpeno_El_Aguacatito_Mex_:250358073
## 3 13FL990008:250304338 171_10_Tuxpeno_El_Aguacatito_Mex_:250358102
## 4 13FL990009:250304339  170_4_Tuxpeno_El_Aguacatito_Mex_:250358089
## 5 13FL990010:250304340  167_5_Tuxpeno_El_Aguacatito_Mex_:250358073
## 6 13FL990012:250304343  166_9_Tuxpeno_El_Aguacatito_Mex_:250358067
##                                        Father     Type
## 1  170_1_Tuxpeno_El_Aguacatito_Mex_:250358087 outcross
## 2 164_14_Tuxpeno_El_Aguacatito_Mex_:250426648 outcross
## 3  164_8_Tuxpeno_El_Aguacatito_Mex_:250358054 outcross
## 4  171_5_Tuxpeno_El_Aguacatito_Mex_:250358097 outcross
## 5  164_1_Tuxpeno_El_Aguacatito_Mex_:250426647 outcross
## 6  164_1_Tuxpeno_El_Aguacatito_Mex_:250426647 outcross
landraces.xo<-read.table(file="xo_Combined_LR13_LR14_parents2_AGPv4_filter0219.txt",
               header=TRUE,stringsAsFactors=FALSE)
head(landraces.xo)
##                  taxon chr   start     end
## 1 13FL990708:250295003   1  262738 1023028
## 2 13FL991485:250349616   1  262855 1030583
## 3 14FL991631:250410929   1  262855 1218925
## 4 14FL991263:250410523   1 1030583 1228056
## 5 14FL990731:250409914   1  262855 1508268
## 6 14FL991166:250410414   1 1228478 1508268

Transforming maize landraces data

#landraces_xo_by_taxon<-landraces.xo %>% 
#  group_by(taxon) %>% 
 #   summarise(xo_by_taxon_count = length(taxon))

#write.table(landraces_xo_by_taxon,file="landraces_xo_by_taxon.txt",sep="\t",quote = #FALSE,row.names = FALSE)


landraces.merged.xo<-merge(landraces.xo,landraces.parentage[1:3],by.x="taxon",by.y="plant_ID")

landraces.child.counts<- landraces.parentage %>% 
  group_by(Mother) %>%
  summarise(family_progeny_count = length(Mother))

landraces.selfed.counts<-landraces.parentage[landraces.parentage$Mother==landraces.parentage$Father,] %>%
   group_by(Mother) %>%
  summarise(family_selfed_count = length(Mother))


landraces.xo.by.mom<- landraces.merged.xo %>% 
  group_by(Mother) %>%
  summarise(family_xo_count = length(Mother))

landraces.merged.xo<-inner_join(landraces.merged.xo,landraces.child.counts)
landraces.merged.xo<-inner_join(landraces.merged.xo,landraces.selfed.counts)
landraces.merged.xo<-inner_join(landraces.merged.xo,landraces.xo.by.mom)


landraces.merged.moms<-inner_join(landraces.xo.by.mom,landraces.child.counts)
landraces.merged.moms<-inner_join(landraces.merged.moms,landraces.selfed.counts)

head(landraces.merged.moms)
## # A tibble: 6 × 4
##                                        Mother family_xo_count
##                                         <chr>           <int>
## 1  164_1_Tuxpeno_El_Aguacatito_Mex_:250426647            4931
## 2 164_11_Tuxpeno_El_Aguacatito_Mex_:250358056            4167
## 3 164_14_Tuxpeno_El_Aguacatito_Mex_:250426648            4972
## 4  164_2_Tuxpeno_El_Aguacatito_Mex_:250358051            5353
## 5  164_4_Tuxpeno_El_Aguacatito_Mex_:250358052            4995
## 6  164_8_Tuxpeno_El_Aguacatito_Mex_:250358054            4827
## # ... with 2 more variables: family_progeny_count <int>,
## #   family_selfed_count <int>
dim(landraces.merged.moms)
## [1] 45  4
landraces.xo.over.progeny.count<-data.frame(rep("landraces",dim(landraces.merged.moms)[1]),landraces.merged.moms$family_xo_count/landraces.merged.moms$family_progeny_count)
names(landraces.xo.over.progeny.count)<-c("teo.or.maize","xo_over_family_count")

head(landraces.xo.over.progeny.count)
##   teo.or.maize xo_over_family_count
## 1    landraces             28.17714
## 2    landraces             30.86667
## 3    landraces             25.49744
## 4    landraces             27.88021
## 5    landraces             28.54286
## 6    landraces             28.39412

About to plot

dim(teosinte.xo.over.progeny.count)
## [1] 52  2
dim(landraces.xo.over.progeny.count)
## [1] 45  2
xo.over.progeny.count<-rbind(teosinte.xo.over.progeny.count,landraces.xo.over.progeny.count)
head(xo.over.progeny.count)
##   teo.or.maize xo_over_family_count
## 1     teosinte             26.58095
## 2     teosinte             23.90667
## 3     teosinte             26.67857
## 4     teosinte             26.20930
## 5     teosinte             26.42391
## 6     teosinte             32.41096
tail(xo.over.progeny.count)
##    teo.or.maize xo_over_family_count
## 92    landraces              6.00000
## 93    landraces             30.21164
## 94    landraces             28.36224
## 95    landraces             28.60000
## 96    landraces             28.03046
## 97    landraces             31.03030
boxplot(xo.over.progeny.count$xo_over_family_count ~xo.over.progeny.count$teo.or.maize,
        main="Per family xo/progeny.count")

Filtering for family.selfed.count>20

  teosinte.merged.moms<-teosinte.merged.moms[teosinte.merged.moms$family_selfed_count>20,]
  teosinte.xo.over.progeny.count<-data.frame(rep("teosinte",dim(teosinte.merged.moms)[1]),teosinte.merged.moms$family_xo_count/teosinte.merged.moms$family_progeny_count)
names(teosinte.xo.over.progeny.count)<-c("teo.or.maize","xo_over_family_count")

  landraces.merged.moms<-landraces.merged.moms[landraces.merged.moms$family_selfed_count>20,]
landraces.xo.over.progeny.count<-data.frame(rep("landraces",dim(landraces.merged.moms)[1]),landraces.merged.moms$family_xo_count/landraces.merged.moms$family_progeny_count)
names(landraces.xo.over.progeny.count)<-c("teo.or.maize","xo_over_family_count")

dim(teosinte.xo.over.progeny.count)
## [1] 42  2
dim(landraces.xo.over.progeny.count)
## [1] 24  2

And plot again

xo.over.progeny.count<-rbind(teosinte.xo.over.progeny.count,landraces.xo.over.progeny.count)
head(xo.over.progeny.count)
##   teo.or.maize xo_over_family_count
## 1     teosinte             26.58095
## 2     teosinte             23.90667
## 3     teosinte             26.67857
## 4     teosinte             26.20930
## 5     teosinte             32.41096
## 6     teosinte             25.35714
tail(xo.over.progeny.count)
##    teo.or.maize xo_over_family_count
## 61    landraces             26.93103
## 62    landraces             29.24194
## 63    landraces             30.21164
## 64    landraces             28.36224
## 65    landraces             28.60000
## 66    landraces             28.03046
boxplot(xo.over.progeny.count$xo_over_family_count ~xo.over.progeny.count$teo.or.maize,
        main="Per family xo/progeny.count")

Getting the medians

landraces.values<-xo.over.progeny.count[xo.over.progeny.count$teo.or.maize=="landraces",]$xo_over_family_count
median(landraces.values)
## [1] 28.26969
teosinte.values<-xo.over.progeny.count[xo.over.progeny.count$teo.or.maize=="teosinte",]$xo_over_family_count
median(teosinte.values)
## [1] 26.53205
wilcox.test(landraces.values,teosinte.values)
## 
##  Wilcoxon rank sum test
## 
## data:  landraces.values and teosinte.values
## W = 705, p-value = 0.006877
## alternative hypothesis: true location shift is not equal to 0

Filtering xo events, removing families with <20 selfed progeny and of the three landrace families with outlier xo/progeny.count

dim(teosinte.merged.xo)
## [1] 127546     10
dim(landraces.merged.xo)
## [1] 137100      9
#number of teosinte families
length(unique(teosinte.merged.xo$Mother))
## [1] 52
#number of landrace families
length(unique(landraces.merged.xo$Mother))
## [1] 45
teosinte.merged.xo<-teosinte.merged.xo[teosinte.merged.xo$family_selfed_count>20,]
landraces.merged.xo<-landraces.merged.xo[landraces.merged.xo$family_selfed_count>20,]


dim(teosinte.merged.xo)
## [1] 109226     10
dim(landraces.merged.xo)
## [1] 101542      9
#number of teosinte families left after filtering
length(unique(teosinte.merged.xo$Mother))
## [1] 42
#number of landrace families left after filtering
length(unique(landraces.merged.xo$Mother))
## [1] 24

Making Maps!!

working with all 10 teosinte chromosomes (20 Mb resolution)

 xo<-teosinte.merged.xo
 num.individuals<-length(unique(xo$taxon)) 
 num.individuals
## [1] 3999
 window.size<-2e7 
 chromo.keep<-c()
 position.keep<-c()
 cM.keep<-c() 
 for (current.chr in 1:10) {

   max.end<-max(filter(xo,chr==current.chr)$end)
   
   my.breaks<-seq(from=0,to=max.end,by=window.size)
   if(max(my.breaks)<max.end){ my.breaks<-c(my.breaks,max.end+1)} #added +1 to max.end to avoid losing last recomb event
    current.chr.table<-hist(xo[xo$chr==current.chr,]$end,breaks = my.breaks,plot=FALSE)  
   
   prob.recombination<- 
   current.chr.table$counts/(2*num.individuals)
   
   #pos<- 50 * log(1/(1-2*prob.recombination))
   pos<-prob.recombination*100
   for (i in 2:length(pos)) {
    pos[i]=pos[i]+pos[i-1]
   }
   pos<-c(0,pos) #adding an start of the chromosome
   chromo.keep<-c(chromo.keep,rep(current.chr,length(pos)))
   cM.keep<-c(cM.keep,pos)
   position.keep<-c(position.keep,my.breaks)
   
 }
 
  roughTeosinteData<-data.frame(chr=chromo.keep,pos=cM.keep,phys_pos=position.keep)
  
  rownames(roughTeosinteData) <- paste0(as.character(roughTeosinteData$chr),"_",round(roughTeosinteData$phys_pos/1e6,digits=0),"Mb")
  head(roughTeosinteData)
##         chr      pos phys_pos
## 1_0Mb     1  0.00000    0e+00
## 1_20Mb    1 37.83446    2e+07
## 1_40Mb    1 53.92598    4e+07
## 1_60Mb    1 65.60390    6e+07
## 1_80Mb    1 74.51863    8e+07
## 1_100Mb   1 80.03251    1e+08
  tail(roughTeosinteData)
##          chr       pos  phys_pos
## 10_60Mb   10  34.22106  60000000
## 10_80Mb   10  37.65941  80000000
## 10_100Mb  10  42.06052 100000000
## 10_120Mb  10  45.37384 120000000
## 10_140Mb  10  64.71618 140000000
## 10_150Mb  10 101.22531 150357240

And now for the landraces (20 Mb resolution)

 xo<-landraces.merged.xo
 num.individuals<-length(unique(xo$taxon)) 
 num.individuals
## [1] 3578
  chromo.keep<-c()
 position.keep<-c()
 cM.keep<-c() 
 for (current.chr in 1:10) {

   max.end<-max(filter(xo,chr==current.chr)$end)
   
   my.breaks<-seq(from=0,to=max.end,by=window.size)
   if(max(my.breaks)<max.end){ my.breaks<-c(my.breaks,max.end)}
   current.chr.table<-hist(xo[xo$chr==current.chr,]$end,breaks = my.breaks,plot=FALSE)  
   
   prob.recombination<- 
   current.chr.table$counts/(2*num.individuals)
   
   #pos<- 50 * log(1/(1-2*prob.recombination))
   pos<-prob.recombination*100
   for (i in 2:length(pos)) {
    pos[i]=pos[i]+pos[i-1]
   }
   pos<-c(0,pos) #adding an start of the chromosome
   chromo.keep<-c(chromo.keep,rep(current.chr,length(pos)))
   cM.keep<-c(cM.keep,pos)
   position.keep<-c(position.keep,my.breaks)
   
 }
 
  roughMaizeLandRaceData<-data.frame(chr=chromo.keep,pos=cM.keep,phys_pos=position.keep)
  
  rownames(roughMaizeLandRaceData) <- paste0(as.character(roughMaizeLandRaceData$chr),"_",round(roughMaizeLandRaceData$phys_pos/1e6,digits=0),"Mb")
  head(roughMaizeLandRaceData)
##         chr      pos phys_pos
## 1_0Mb     1  0.00000    0e+00
## 1_20Mb    1 36.76635    2e+07
## 1_40Mb    1 54.24818    4e+07
## 1_60Mb    1 69.07490    6e+07
## 1_80Mb    1 76.46730    8e+07
## 1_100Mb   1 82.85355    1e+08
  tail(roughMaizeLandRaceData)
##          chr       pos  phys_pos
## 10_60Mb   10  47.47065  60000000
## 10_80Mb   10  50.29346  80000000
## 10_100Mb  10  54.38793 100000000
## 10_120Mb  10  57.79765 120000000
## 10_140Mb  10  77.36165 140000000
## 10_151Mb  10 115.42761 150827693

Total map length in cM

  ##for teosinte
  max.per.chromosome.teo<-tapply(roughTeosinteData$pos,roughTeosinteData$chr,max)
  max.per.chromosome.teo
##        1        2        3        4        5        6        7        8 
## 200.7627 166.9542 141.1353 132.7332 143.2608 110.3651 131.1953 126.0065 
##        9       10 
## 112.0280 101.2253
  total.cM.teosinte<-sum(max.per.chromosome.teo)
  total.cM.teosinte
## [1] 1365.666
  ## for landraces
  max.per.chromosome.landraces<-tapply(roughMaizeLandRaceData$pos,roughMaizeLandRaceData$chr,max)
  max.per.chromosome.landraces
##        1        2        3        4        5        6        7        8 
## 192.1465 162.5349 147.6942 144.7736 148.4069 117.5377 147.4846 125.3913 
##        9       10 
## 117.5797 115.4276
  total.cM.landraces<-sum(max.per.chromosome.landraces)
  total.cM.landraces
## [1] 1418.977

Chromosome by chromosome comparisons

teosinteMap<- table2map(roughTeosinteData)
landracesMap<- table2map(roughMaizeLandRaceData)
for(my.chr in 1:10) {
  par(mfrow=c(1,2))
  plotMap(teosinteMap,chr=c(my.chr),show.marker.names=TRUE,main=paste0("Chromosome ",my.chr ," teosinte"))
  plotMap(landracesMap,chr=c(my.chr), show.marker.names=TRUE, main=paste0("Chromosome ",my.chr ,"  landraces"))
  
  
}

Making the maps again but at 1Mb resolution

working with all 10 teosinte chromosomes (1 Mb resolution)

 xo<-teosinte.merged.xo
 num.individuals<-length(unique(xo$taxon)) 
 num.individuals
## [1] 3999
 window.size<-1e6  #1Mb 
 chromo.keep<-c()
 position.keep<-c()
 cM.keep<-c() 
 for (current.chr in 1:10) {

   max.end<-max(filter(xo,chr==current.chr)$end)
   
   my.breaks<-seq(from=0,to=max.end,by=window.size)
   if(max(my.breaks)<max.end){ my.breaks<-c(my.breaks,max.end)}
   current.chr.table<-hist(xo[xo$chr==current.chr,]$end,breaks = my.breaks,plot=FALSE)  
   
   prob.recombination<- 
   current.chr.table$counts/(2*num.individuals)
   
   #pos<- 50 * log(1/(1-2*prob.recombination))
   pos<-prob.recombination*100
   for (i in 2:length(pos)) {
    pos[i]=pos[i]+pos[i-1]
   }
   pos<-c(0,pos) #adding an start of the chromosome
   chromo.keep<-c(chromo.keep,rep(current.chr,length(pos)))
   cM.keep<-c(cM.keep,pos)
   position.keep<-c(position.keep,my.breaks)
   
 }
 
  fineTeosinteData<-data.frame(chr=chromo.keep,pos=cM.keep,phys_pos=position.keep)
  
  rownames(fineTeosinteData) <- paste0(as.character(fineTeosinteData$chr),"_",round(fineTeosinteData$phys_pos/1e6,digits=2),"Mb")
  head(fineTeosinteData)
##       chr        pos phys_pos
## 1_0Mb   1  0.0000000    0e+00
## 1_1Mb   1  0.0500125    1e+06
## 1_2Mb   1  0.6626657    2e+06
## 1_3Mb   1  2.5006252    3e+06
## 1_4Mb   1  7.2518130    4e+06
## 1_5Mb   1 10.2775694    5e+06
  tail(fineTeosinteData)
##             chr       pos  phys_pos
## 10_146Mb     10  83.25831 146000000
## 10_147Mb     10  87.18430 147000000
## 10_148Mb     10  90.72268 148000000
## 10_149Mb     10  96.91173 149000000
## 10_150Mb     10 101.12528 150000000
## 10_150.36Mb  10 101.22531 150357239

And now for the landraces (1 Mb resolution)

 xo<-landraces.merged.xo
 num.individuals<-length(unique(xo$taxon)) 
 num.individuals
## [1] 3578
  chromo.keep<-c()
 position.keep<-c()
 cM.keep<-c() 
 for (current.chr in 1:10) {

   max.end<-max(filter(xo,chr==current.chr)$end)
   
   my.breaks<-seq(from=0,to=max.end,by=window.size)
   if(max(my.breaks)<max.end){ my.breaks<-c(my.breaks,max.end)}
   current.chr.table<-hist(xo[xo$chr==current.chr,]$end,breaks = my.breaks,plot=FALSE)  
   
   prob.recombination<- 
   current.chr.table$counts/(2*num.individuals)
   
   #pos<- 50 * log(1/(1-2*prob.recombination))
   pos<-prob.recombination*100
   for (i in 2:length(pos)) {
    pos[i]=pos[i]+pos[i-1]
   }
   pos<-c(0,pos) #adding an start of the chromosome
   chromo.keep<-c(chromo.keep,rep(current.chr,length(pos)))
   cM.keep<-c(cM.keep,pos)
   position.keep<-c(position.keep,my.breaks)
   
 }
 
  fineMaizeLandRaceData<-data.frame(chr=chromo.keep,pos=cM.keep,phys_pos=position.keep)
  
  rownames(fineMaizeLandRaceData) <- paste0(as.character(fineMaizeLandRaceData$chr),"_",round(fineMaizeLandRaceData$phys_pos/1e6,digits=2),"Mb")
  head(fineMaizeLandRaceData)
##       chr       pos phys_pos
## 1_0Mb   1 0.0000000    0e+00
## 1_1Mb   1 0.0000000    1e+06
## 1_2Mb   1 0.5310229    2e+06
## 1_3Mb   1 2.2498603    3e+06
## 1_4Mb   1 6.3722750    4e+06
## 1_5Mb   1 9.2230296    5e+06
  tail(fineMaizeLandRaceData)
##             chr       pos  phys_pos
## 10_146Mb     10  97.59642 146000000
## 10_147Mb     10 100.72666 147000000
## 10_148Mb     10 105.88317 148000000
## 10_149Mb     10 111.15148 149000000
## 10_150Mb     10 115.31582 150000000
## 10_150.83Mb  10 115.42761 150827693

Total map length in cM

  ##for teosinte
  max.per.chromosome.teo<-tapply(fineTeosinteData$pos,fineTeosinteData$chr,max)
  max.per.chromosome.teo
##        1        2        3        4        5        6        7        8 
## 200.7627 166.9542 141.1353 132.7332 143.2608 110.3651 131.1953 126.0065 
##        9       10 
## 112.0280 101.2253
  total.cM.teosinte<-sum(max.per.chromosome.teo)
  total.cM.teosinte
## [1] 1365.666
  ## for landraces
  max.per.chromosome.landraces<-tapply(fineMaizeLandRaceData$pos,fineMaizeLandRaceData$chr,max)
  max.per.chromosome.landraces
##        1        2        3        4        5        6        7        8 
## 192.1465 162.5349 147.6942 144.7736 148.4069 117.5377 147.4846 125.3913 
##        9       10 
## 117.5797 115.4276
  total.cM.landraces<-sum(max.per.chromosome.landraces)
  total.cM.landraces
## [1] 1418.977

Chromosome by chromosome comparisons

fineTeosinteMap<- table2map(fineTeosinteData)
fineLandracesMap<- table2map(fineMaizeLandRaceData)
for(my.chr in 1:10) {
  par(mfrow=c(1,2))
  plotMap(fineTeosinteMap,chr=c(my.chr),show.marker.names=TRUE,main=paste0("Chromosome ",my.chr ," teosinte"))
  plotMap(fineLandracesMap,chr=c(my.chr), show.marker.names=TRUE, main=paste0("Chromosome ",my.chr ,"  landraces"))
  
  
}

## Formatting 1MB Maize landraces map for eagle2

  head(fineMaizeLandRaceData)
##       chr       pos phys_pos
## 1_0Mb   1 0.0000000    0e+00
## 1_1Mb   1 0.0000000    1e+06
## 1_2Mb   1 0.5310229    2e+06
## 1_3Mb   1 2.2498603    3e+06
## 1_4Mb   1 6.3722750    4e+06
## 1_5Mb   1 9.2230296    5e+06
  names(fineMaizeLandRaceData)
## [1] "chr"      "pos"      "phys_pos"
  names(fineMaizeLandRaceData)[2]<-"Genetic_Map(cM)"
  names(fineMaizeLandRaceData)[3]<-"pos"
  head(fineMaizeLandRaceData)
##       chr Genetic_Map(cM)   pos
## 1_0Mb   1       0.0000000 0e+00
## 1_1Mb   1       0.0000000 1e+06
## 1_2Mb   1       0.5310229 2e+06
## 1_3Mb   1       2.2498603 3e+06
## 1_4Mb   1       6.3722750 4e+06
## 1_5Mb   1       9.2230296 5e+06
  combined<-c()
  for(i in 1:(dim(fineMaizeLandRaceData)[1]-1)){
    if(!(fineMaizeLandRaceData$chr[i]==fineMaizeLandRaceData$chr[i+1])) {
      combined[i]<-0
    
    } else {
     combined[i]<-(fineMaizeLandRaceData$`Genetic_Map(cM)`[i+1]-fineMaizeLandRaceData$`Genetic_Map(cM)`[i]) /
                  ((fineMaizeLandRaceData$pos[i+1]-fineMaizeLandRaceData$pos[i])/1e6)   
    }
  }

 combined[dim(fineMaizeLandRaceData)[1]]<-0
 fineMaizeLandRaceData.4.Eagle2<-data.frame(fineMaizeLandRaceData$chr,
                                       fineMaizeLandRaceData$pos,
                                       combined,
                                       fineMaizeLandRaceData$`Genetic_Map(cM)`)
 
  names(fineMaizeLandRaceData.4.Eagle2)<-c("chr","position","COMBINED_rate(cM/Mb)","Genetic_Map(cM)")    
  head(fineMaizeLandRaceData.4.Eagle2) 
##   chr position COMBINED_rate(cM/Mb) Genetic_Map(cM)
## 1   1    0e+00            0.0000000       0.0000000
## 2   1    1e+06            0.5310229       0.0000000
## 3   1    2e+06            1.7188373       0.5310229
## 4   1    3e+06            4.1224148       2.2498603
## 5   1    4e+06            2.8507546       6.3722750
## 6   1    5e+06            5.1145892       9.2230296
  write.table(format(fineMaizeLandRaceData.4.Eagle2,digits=5),file="landraces_mapForEagle2_1Mb.txt",row.names = FALSE,quote = FALSE )

Formatting 1MB teosinte map for eagle2

  head(fineTeosinteData)
##       chr        pos phys_pos
## 1_0Mb   1  0.0000000    0e+00
## 1_1Mb   1  0.0500125    1e+06
## 1_2Mb   1  0.6626657    2e+06
## 1_3Mb   1  2.5006252    3e+06
## 1_4Mb   1  7.2518130    4e+06
## 1_5Mb   1 10.2775694    5e+06
  names(fineTeosinteData)
## [1] "chr"      "pos"      "phys_pos"
  names(fineTeosinteData)[2]<-"Genetic_Map(cM)"
  names(fineTeosinteData)[3]<-"pos"
  head(fineTeosinteData)
##       chr Genetic_Map(cM)   pos
## 1_0Mb   1       0.0000000 0e+00
## 1_1Mb   1       0.0500125 1e+06
## 1_2Mb   1       0.6626657 2e+06
## 1_3Mb   1       2.5006252 3e+06
## 1_4Mb   1       7.2518130 4e+06
## 1_5Mb   1      10.2775694 5e+06
  combined<-c()
  for(i in 1:(dim(fineTeosinteData)[1]-1)){
    if(!(fineTeosinteData$chr[i]==fineTeosinteData$chr[i+1])) {
      combined[i]<-0
    
    } else {
     combined[i]<-(fineTeosinteData$`Genetic_Map(cM)`[i+1]-fineTeosinteData$`Genetic_Map(cM)`[i]) /
                  ((fineTeosinteData$pos[i+1]-fineTeosinteData$pos[i])/1e6)   
    }
  }

 combined[dim(fineTeosinteData)[1]]<-0
 fineTeosinteData.4.Eagle2<-data.frame(fineTeosinteData$chr,
                                       fineTeosinteData$pos,
                                       combined,
                                       fineTeosinteData$`Genetic_Map(cM)`)
 
  names(fineTeosinteData.4.Eagle2)<-c("chr","position","COMBINED_rate(cM/Mb)","Genetic_Map(cM)")    
  head(fineTeosinteData.4.Eagle2) 
##   chr position COMBINED_rate(cM/Mb) Genetic_Map(cM)
## 1   1    0e+00            0.0500125       0.0000000
## 2   1    1e+06            0.6126532       0.0500125
## 3   1    2e+06            1.8379595       0.6626657
## 4   1    3e+06            4.7511878       2.5006252
## 5   1    4e+06            3.0257564       7.2518130
## 6   1    5e+06            5.2388097      10.2775694
  write.table(format(fineTeosinteData.4.Eagle2,digits=5),file="teosinte_mapForEagle2_1Mb.txt",row.names = FALSE,quote = FALSE )

Transforming data for recombination rate plots

##putting eagle2 formatted maps together (teosinte and maize)
 fineMaizeLandRaceData.4.Eagle2<-cbind(fineMaizeLandRaceData.4.Eagle2,species=rep("maize.landrace",dim(fineMaizeLandRaceData.4.Eagle2)[1]))
head(fineMaizeLandRaceData.4.Eagle2)
##   chr position COMBINED_rate(cM/Mb) Genetic_Map(cM)        species
## 1   1    0e+00            0.0000000       0.0000000 maize.landrace
## 2   1    1e+06            0.5310229       0.0000000 maize.landrace
## 3   1    2e+06            1.7188373       0.5310229 maize.landrace
## 4   1    3e+06            4.1224148       2.2498603 maize.landrace
## 5   1    4e+06            2.8507546       6.3722750 maize.landrace
## 6   1    5e+06            5.1145892       9.2230296 maize.landrace
fineTeosinteData.4.Eagle2<-cbind(fineTeosinteData.4.Eagle2,species=rep("teosinte",dim(fineTeosinteData.4.Eagle2)[1]))
head(fineTeosinteData.4.Eagle2)
##   chr position COMBINED_rate(cM/Mb) Genetic_Map(cM)  species
## 1   1    0e+00            0.0500125       0.0000000 teosinte
## 2   1    1e+06            0.6126532       0.0500125 teosinte
## 3   1    2e+06            1.8379595       0.6626657 teosinte
## 4   1    3e+06            4.7511878       2.5006252 teosinte
## 5   1    4e+06            3.0257564       7.2518130 teosinte
## 6   1    5e+06            5.2388097      10.2775694 teosinte
teo.and.landraces<-rbind(fineMaizeLandRaceData.4.Eagle2,fineTeosinteData.4.Eagle2)

for(chr.index in 1:10){

 print(
   ggplot(data=filter(teo.and.landraces,chr==chr.index),aes(x=position,y=`COMBINED_rate(cM/Mb)`))+geom_line(mapping=aes(color=species)) +labs(title=paste0("Chr ",chr.index))
        
 )
}

for(chr.index in 1:10){

 print(
   ggplot(data=filter(teo.and.landraces,chr==chr.index),aes(x=position,y=`Genetic_Map(cM)`))+geom_line(mapping=aes(color=species)) +labs(title=paste0("Chr ",chr.index))
        
 )
}