# This skript is for the analysis of QTL - Quantitative Trait Locus
#
# get access to hyper data; look at its class
rm(list=ls())
library(qtl) # a special package created for qtl analysis
data(hyper)
###### Let's now explore the dataset
summary(hyper) # gives a summary
help(hyper)
names(hyper)
# let's have a deeper look into the object genotype
names(hyper$geno)
names(hyper$geno$`1`)
hyper$geno$`1`$data
hyper$geno$`1`$map
# let's have a deeper look into genotype
names(hyper$pheno)
hyper$pheno
hyper$pheno[1:5,]
summary(hyper$pheno$bp)
## some more interesting things to ask to the dataset with built-in functions:
nphe(hyper) # Number of phenotypes
nind(hyper) # number of individuals
nchr(hyper) # number of Chromosomes
# now have a more graphical look on the dataset
plot.map(hyper)
# if you may want to save this graphic for later, you could easily store
# plots somewhere with the following commands:
jpeg("geneticmap.jpg") # this tells R to create a new imagefile
plot.map(hyper)
graphics.off() # this tells R to close the imagefile and to store it in the working directory
getwd() # if you don't know your working directory, this command will tell you
# and compare it with the data we have.
#
# A diagnostic tool
# calculationg there any genotypes that accociated with an statistical error
# with the probability of 0.01
hyper <-calc.errorlod(hyper, error.prob=0.01)
top.errorlod(hyper)
hyper <- est.rf(hyper) # creates a recombination fraction chart
hyper
plot.geno(hyper,chr=3, 30:150)
newmap <- est.map(hyper) ## creats a so-called "estimated map".
#The idea behind the "estimated map" is:
plot.map(hyper,newmap) # Based on the data we have there le'ts create a genetic map that we would EXPECT to see
plot.rf(hyper)
help(plot.geno)
hyper <-clean(hyper) ## cleans calculations done on t
# the next step (after inspecting the data) is to calculate the genotype-probability
# The command "calc.genoprob" is using the HMM-Algorithmus of Haldane for a probability analysis of genotypes
# in your dataset. And simulates a dataset and various possible permutations of markers for statistical tests
GT_prob <-calc.genoprob(hyper,step=2,off.end = 0.0,error.prob=1.0e-4,map.function="haldane",stepwidth = "fixed")
hyper <- sim.geno(hyper,n.draws=32,step=2,off.end = 0.0,error.prob = 1.0e-4,map.function="haldane",stepwidth = "fixed")
# Scan for associations: the command "scanone" performs a genome scan with a single QTL model,
# with possible allowance for covariates, using any of several possible models for the phenotype
# and any of several possible numerical methods.
hyper.scan <- scanone(hyper, pheno.col = 1, model="normal",method="em")
hyper.scan.perm <- scanone(hyper, pheno.col = 1, model="normal")
# define threshold-values in order to compare our data
thresh <- summary(hyper.scan.perm)
head(thresh)
# Now we are ready to have a look at predicted QTL
#
plot(hyper.scan,main="Plot of possible QTL")
summary(hyper.scan)
#Let's find and analyse now the interesting marker in our dataset and plot them against the
# Gorresponding allel
Marker_name1 <- find.marker(hyper,chr=1,pos=49.2)
effectplot(hyper, pheno.col = 1,mname1=Marker_name1 )
Marker_name2 <- find.marker(hyper,chr=4,pos=51.7)
effectplot(hyper, pheno.col = 1,mname1=Marker_name2 )
# Now let's calculate a confidence interval whether or not we can have confidence in our data
Confidence_Chr1 <- bayesint(hyper.scan,chr = 1, prob = 0.95)
Confidence_Chr1
plot(hyper.scan,chr = 1,lodcolumn = 1,main="Confidence interval for Chromosom Nr. 1")
lines(x=Confidence_Chr1[c(1,3),2], y=c(0,0),type ="l",col="green",lwd=5)
# let's do the same for Chromosome Nr. 2
Confidence_Chr2 <- bayesint(hyper.scan,chr = 2, prob = 0.95)
Confidence_Chr2
plot(hyper.scan,chr = 2,lodcolumn = 1,main="Confidence interval for Chromosom Nr. 2")
lines(x=Confidence_Chr2[c(1,3),2], y=c(0,0),type ="l",col="green",lwd=5)
# let's do the same for Chromosome Nr. 4
Confidence_Chr4 <- bayesint(hyper.scan,chr = 4, prob = 0.95)
Confidence_Chr4
plot(hyper.scan,chr = 4,lodcolumn = 1,main="Confidence interval for Chromosom Nr. 4")
lines(x=Confidence_Chr4[c(1,3),2], y=c(0,0),type ="l",col="green",lwd=5)
# let's do the same for Chromosome Nr. 15
Confidence_Chr15 <- bayesint(hyper.scan,chr = 15, prob = 0.95)
Confidence_Chr15
plot(hyper.scan,chr = 15,lodcolumn = 1,main="Confidence interval for Chromosom Nr. 15")
lines(x=Confidence_Chr15[c(1,3),2], y=c(0,0),type ="l",col="green",lwd=5)
plot(hyper.scan,main="Result of scan for QTL")
jpeg("QTL-Regions.jpg") # this tells R to create a new imagefile
plot(hyper.scan,main="Result of scan for QTL")
graphics.off() # this tells R to close the imagefile and to store it in the working directory