Skip to content
Snippets Groups Projects
Commit ea06b0f1 authored by Marie Bourdon's avatar Marie Bourdon
Browse files

mark_prop chi2 N2

parent 612b2bc8
No related branches found
No related tags found
No related merge requests found
place <- c()
previous <- c()
follow <- c()
#get information in newmap
for(i in names(map)){
marker <- c(marker,names(map[[i]]))
chr <- c(chr,rep(i,times=length(map[[i]])))
pos <- c(pos,unname(map[[i]]))
place <- c(place,"first",rep("middle",times=(length(map[[i]])-2)),"last")
prev <- c(NA,unname(map[[i]])[1:length(map[[i]])-1])
previous <- c(previous,prev)
fol <- c(unname(map[[i]])[2:length(map[[i]])],NA)
follow <- c(follow,fol)
}
tab_map <- tibble(marker = marker,
chr = chr,
place = place,
pos = pos,
previous = pos-previous,
follow = follow-pos)
tab_map <- tab_map %>% mutate(exclude=case_when(place == "first" & follow > 100 ~ 1,
place == "middle" & previous > 100 & follow > 100 ~ 1,
place == "last" & previous > 100 ~ 1,
T ~ 0))
bad_markers <- tab_map %>% filter(exclude == 1) %>% pull(marker) bad_markers <- tab_map %>% filter(exclude == 1) %>% pull(marker)
} }
mark_estmap <- function(map,dist=100){ mark_estmap <- function(map,dist=100){
...@@ -510,3 +486,27 @@ rm(stuart_cross) ...@@ -510,3 +486,27 @@ rm(stuart_cross)
library(qtl) library(qtl)
data(stuart_cross) data(stuart_cross)
summary(stuart_cross) summary(stuart_cross)
setwd("~/Documents/stuart_package/stuart")
devtools::build(path=".",vignettes=FALSE)
devtools::build_vignettes()
devtools::build_manual()
library(stuart)
mark_prop
devtools::build_manual(path=".")
library(dplyr)
library(stuart)
annot_mini <- read.csv(url("https://raw.githubusercontent.com/kbroman/MUGAarrays/master/UWisc/mini_uwisc_v2.csv"))
data(genos)
summary(genos)
data(phenos)
summary(phenos)
strains <- geno_strains(ref=annot_mini,geno=genos,par1=c("StrainsA_1","StrainsA_2"),
par2=c("StrainsB_1","StrainsB_2"),name1="parent1",name2="parent2")
head(strains) %>% print.data.frame()
genos <- genos %>% filter(!Sample.ID %in% c("StrainsA_1", "StrainsA_2",
"StrainsB_1","StrainsB_2"))
# how to use the function:
# stuart_tab <- tab_mark(genos)
data(stuart_tab)
summary(stuart_tab)
View(stuart_tab)
...@@ -61,27 +61,36 @@ mark_prop <- function(tab,cross,homo=NA,hetero=NA,pval=NA,na=0.5){ ...@@ -61,27 +61,36 @@ mark_prop <- function(tab,cross,homo=NA,hetero=NA,pval=NA,na=0.5){
T ~ exclude_prop T ~ exclude_prop
)) ))
#exclude with pval chisq.test #exclude with pval chisq.test
} else if(is.na(pval)==FALSE){ } else if(is.na(pval)==FALSE){
#if cross F2 #if cross F2
if(cross=="F2"){ if(cross=="F2"){
#chisq test
tab <- tab %>% filter(p_NA != 1) %>% rowwise() %>% tab <- tab %>% filter(p_NA != 1) %>% rowwise() %>%
mutate(.,chi_pval = tibble(n_HM1,n_HM2,n_HT) %>% mutate(.,chi_pval = tibble(n_HM1,n_HM2,n_HT) %>%
chisq.test(p=c(0.25,0.25,0.5)) %>% .$p.value) %>% chisq.test(p=c(0.25,0.25,0.5)) %>% .$p.value) %>%
full_join(.,tab,by=all_of(names)) full_join(.,tab,by=all_of(names))
#convert in score
tab <- tab %>% mutate(exclude_prop=case_when(chi_pval < pval ~ 1, tab <- tab %>% mutate(exclude_prop=case_when(chi_pval < pval ~ 1,
T ~ exclude_prop)) T ~ exclude_prop))
#if cross N2 #if cross N2
} else if(cross=="N2"){ } else if(cross=="N2"){
#stock homozygous data in a new column (either n_HM1 or n-<hm2)
tab <- tab %>% mutate(n_hm = case_when(n_HM1>n_HM2 ~ n_HM1,
T ~ n_HM2))
#chisq test
tab <- tab %>% filter(p_NA != 1) %>% rowwise() %>% tab <- tab %>% filter(p_NA != 1) %>% rowwise() %>%
mutate(.,chi_pval = tibble(n_HM1,n_HT) %>% mutate(.,chi_pval =chisq.test(x=c(n_hm,n_HT),p=c(0.5,0.5)) %>% .$p.value) %>%
chisq.test(p=c(0.5,0.5)) %>% .$p.value) %>% full_join(.,tab,by=c(all_of(names),"n_hm"))
full_join(.,tab,by=all_of(names))
#convert in score
tab <- tab %>% mutate(exclude_prop=case_when(chi_pval < pval ~ 1, tab <- tab %>% mutate(exclude_prop=case_when(chi_pval < pval ~ 1,
(n_HM1 !=0 & n_HM2 != 0) ~ 1,
T ~ exclude_prop)) T ~ exclude_prop))
} }
......
This diff is collapsed.
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment