maps <- readRDS("../2_LinMap/10_4_linkageMap.RData")
library(LinkageMapView)
library(janitor)
##
## Attaching package: 'janitor'
## The following objects are masked from 'package:stats':
##
## chisq.test, fisher.test
library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──
## ✓ ggplot2 3.3.5 ✓ purrr 0.3.4
## ✓ tibble 3.1.6 ✓ dplyr 1.0.7
## ✓ tidyr 1.2.0 ✓ stringr 1.4.0
## ✓ readr 1.4.0 ✓ forcats 0.5.1
## Warning: package 'tidyr' was built under R version 4.0.5
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
Convert to LinkageMapView format
LMV_male <- maps %>%
select(LG,male,SNPID) %>%
drop_na() %>%
rename(group = LG, position = male,
locus = SNPID) %>%
mutate(across(where(is.factor),
as.character)) %>%
mutate(group = paste0(group, "M"))
LMV_female <- maps %>%
select(LG,female,SNPID) %>%
drop_na() %>%
rename(group = LG, position = female,
locus = SNPID) %>%
mutate(across(where(is.factor),
as.character)) %>%
mutate(group = paste0(group, "F"))
LMV_average <- maps %>%
select(LG,sex.averaged,SNPID) %>%
drop_na() %>%
rename(group = LG, position = sex.averaged,
locus = SNPID) %>%
mutate(across(where(is.factor),
as.character)) %>%
mutate(group = paste0(group, "A"))
HC <- rbind(LMV_male, LMV_female, LMV_average)
levels <- c("1F","1M","1A","2F","2M","2A","3F","3M","3A",
"4F","4M","4A","5F","5M","5A","6F","6M","6A",
"7F","7M","7A","8F","8M","8A","9F","9M","9A",
"10F","10M","10A","11F","11M","11A",
"12F","12M","12A","13F","13M","13A")
HC$group <- fct_relevel(HC$group,
levels = levels)
## Warning: Outer names are only allowed for unnamed scalar atomic inputs
HC <- HC[order(HC$group,HC$position),]
outfile = file.path("../2_LinMap/12_1_maleMap_density.pdf")
lmv.linkage.plot(LMV_male,outfile,denmap=TRUE,
pdf.width = 15, pdf.height = 15,
cex.axis = 2,lgw = 0.30,cex.lgtitle = 2,
font.main = 2,cex.main = 3,
main="Male linkage maps")
## Required pdf.width = 11.67
## Required pdf.height = 4.56388888888891
## Using pdf.width = 15
## Using pdf.height = 15
outfile = file.path("../2_LinMap/12_1_femaleMap_density.pdf")
lmv.linkage.plot(LMV_female,outfile,denmap=TRUE,
pdf.width = 15, pdf.height = 15,
cex.axis = 2,lgw = 0.30,cex.lgtitle = 2,
font.main = 2,cex.main = 3,
main="Female linkage maps")
## Required pdf.width = 10.708
## Required pdf.height = 6.23402777777799
## Using pdf.width = 15
## Using pdf.height = 15
outfile = file.path("../2_LinMap/12_1_averageMap_density.pdf")
lmv.linkage.plot(LMV_average,outfile,denmap=TRUE,
pdf.width = 15, pdf.height = 15,
cex.axis = 2,lgw = 0.30,cex.lgtitle = 2,
font.main = 2,cex.main = 3,
main="Sex-averaged linkage maps")
## Required pdf.width = 10.9506666666667
## Required pdf.height = 13.5656250000006
## Using pdf.width = 15
## Using pdf.height = 15
Now let’s plot the male, female and sex-averaged linkage maps showing only the position of markers as black horizontal lines.
outfile = file.path("../2_LinMap/12_2_maleMap.pdf")
lmv.linkage.plot(LMV_male,outfile, showonly = c("22896_41"),
ruler = TRUE,cex.axis = 2,
lgw = 0.30,
lgtitles = c("LG1","LG2","LG3","LG4",
"LG5","LG6","LG7","LG8",
"LG9","LG10","LG11","LG12",
"LG13"),
cex.lgtitle = 2,
pdf.width = 20, pdf.height = 18,
main="Male linkage maps",
font.main = 2,cex.main = 3)
## Required pdf.width = 14.8961666666667
## Required pdf.height = 1.38506666666667
## Using pdf.width = 20
## Using pdf.height = 18
outfile = file.path("../2_LinMap/12_2_femaleMap.pdf")
lmv.linkage.plot(LMV_female,outfile, showonly = c("22896_41"),
ruler = TRUE,cex.axis = 2,
lgw = 0.30,
lgtitles = c("LG1","LG2","LG3","LG4",
"LG5","LG6","LG7","LG8",
"LG9","LG10","LG11","LG12",
"LG13"),
cex.lgtitle = 2,
pdf.width = 20, pdf.height = 18,
main="Female linkage maps",
font.main = 2,cex.main = 3)
## Required pdf.width = 14.8961666666667
## Required pdf.height = 1.5526
## Using pdf.width = 20
## Using pdf.height = 18
outfile = file.path("../2_LinMap/12_2_averageMap.pdf")
lmv.linkage.plot(LMV_average,outfile, showonly = c("22896_41"),
ruler = TRUE,cex.axis = 2,
lgw = 0.30,
lgtitles = c("LG1","LG2","LG3","LG4",
"LG5","LG6","LG7","LG8",
"LG9","LG10","LG11","LG12",
"LG13"),
cex.lgtitle = 2,
pdf.width = 20, pdf.height = 18,
main="Sex-averaged linkage maps",
font.main = 2,cex.main = 3)
## Required pdf.width = 14.8961666666667
## Required pdf.height = 1.5526
## Using pdf.width = 20
## Using pdf.height = 18
First, let’s create a dataframe with male and female maps, and then order them.
temp <- rbind(LMV_male, LMV_female) %>% arrange(group)
temp1 <- temp %>% filter(group != "10M" & group != "11M" &
group != "12M" & group != "13M" &
group != "10F" & group != "11F" &
group != "12F" & group != "13F")
temp2 <- temp %>% filter(group == "10M" | group == "11M" |
group == "12M" | group == "13M" |
group == "10F" | group == "11F" |
group == "12F" | group == "13F")
temp <- rbind(temp1, temp2)
And then plot
outfile = file.path("../2_LinMap/12_3_maleVSfemaleMap.pdf")
lmv.linkage.plot(temp,outfile, showonly = c("22896_41"),
ruler = TRUE,cex.axis = 3,
lgw = 0.30,
lgtitles = c("LG1f","LG1m","LG2f","LG2m",
"LG3f","LG3m","LG4f","LG4m",
"LG5f","LG5m","LG6f","LG6m",
"LG7f","LG7m","LG8f","LG8m",
"LG9f","LG9m","LG10f","LG10m",
"LG11f","LG11m","LG12f","LG12m",
"LG13f","LG13m"),
cex.lgtitle = 3,
pdf.width = 50, pdf.height = 20,
main="Female and male linkage maps",
font.main = 2,cex.main = 4)
## Required pdf.width = 28.6923333333333
## Required pdf.height = 1.7526
## Using pdf.width = 50
## Using pdf.height = 20
And now let’s plot the density map for LG1 for both the female and male next to each other to show the heterochiasmy.
temp <- temp %>% filter(group == "1F" | group == "1M")
outfile = file.path("../2_LinMap/12_4_maleVSfemaleLG1_density.pdf")
lmv.linkage.plot(temp, outfile, denmap=TRUE,
pdf.width = 3, pdf.height = 7,
cex.axis = 1, lgw = 0.30, cex.lgtitle = 1,
font.main = 2, cex.main = 1,
lgtitles = c("Female","Male"),
main="Female and male LG1")
## Required pdf.width = 1.80366666666667
## Required pdf.height = 6.03402777777799
## Using pdf.width = 3
## Using pdf.height = 7
## 5. Plot position in male versus female maps
rm(list = ls())
Read in map data
maps <- readRDS("../2_LinMap/10_4_linkageMap.RData")
Relevel in order of LG name
maps <- maps %>%
mutate(LG = fct_relevel(LG, "1", "2", "3", "4",
"5", "6", "7", "8", "9",
"10", "11", "12", "13"))
And plot
p1 <- ggplot(maps,
aes(x=female, y=male,
group = LG, colour = LG)) +
geom_point(size=0.5) + theme_classic() +
theme(legend.position = "none") +
xlab("Female map position (cM)") +
ylab("Male map position (cM)") +
scale_y_continuous(limits = c(0, 170),
breaks = seq(0, 150, by = 50)) +
scale_x_continuous(limits = c(0,265),
breaks = seq(0, 250, by = 50))
p1
## Warning: Removed 7015 rows containing missing values (geom_point).
Save the resulting plot
ggsave(filename = "../2_LinMap/12_5_MaleVsFemale_All.jpg",
plot = p1, width = 6, height = 4)
## Warning: Removed 7015 rows containing missing values (geom_point).
That figure is a bit confusing, showing many LGs at once. Let’s plot each LG individually.
maps.NA <- maps %>% drop_na()
p2 <- ggplot(maps.NA,
aes(x=female, y=male,
group = LG, colour = LG)) +
geom_point() + facet_wrap(~LG, ncol = 5) +
xlab("Female map position (cM)") +
ylab("Male map position (cM)") +
theme_classic() +
theme(legend.position = "none")
p2
Save the resulting figure
ggsave(filename = "../2_LinMap/12_6_MaleVsFemale.jpg",
plot = p2, width = 8, height = 6)
rm(p1,p2)
Now, let’s plot marker density of male versus female maps to show the high level of heterochiasmy.