Mercurial > repos > richjp > cloudmap_zebrafish_hom_mapping

changeset 9:0fee75e8a18e draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Deleted selected files
author richjp
date Mon, 13 Oct 2014 08:49:26 -0400
parents 0ede96d29bd5
children 9b60d5a2c313
files cloudmap_zf_hom_mapping.xml
diffstat 1 files changed, 0 insertions(+), 1301 deletions(-) [+]
line wrap: on
line diff
--- a/cloudmap_zf_hom_mapping.xml	Mon Oct 13 08:45:03 2014 -0400
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,1301 +0,0 @@
-<tool id="CloudMap: zebrafish homozygosity mapping" name="CloudMap: zebrafish homozygosity mapping" version="0.01">
-<description>Map a zebrafish mutant using bulk sergeant homozygosity mapping</description>
-<requirements>
-      <requirement type="package" version="9.10">ghostscript</requirement>
-      <requirement type="package" version="1.3.18">graphicsmagick</requirement>
-  </requirements>
-<command interpreter="python">
-
-        cloudmap_zf_hom_mapping.py --script_path "$runMe" --interpreter "Rscript" 
-            --tool_name "cloudmap_zf_hom_mapping" --input_tab "$input1"   --output_dir "$html_file.files_path" --output_html "$html_file" --make_HTML "yes" 
-</command>
-<inputs>
-<param name="input1"  type="data" format="vcf" label="Select a suitable input file from your history"/> 
-<param name="job_name" type="text" label="Supply a name for the outputs to remind you what they contain" value="cloudmap_zf_hom_mapping"/> 
-
-</inputs>
-<outputs>
- <data format="html" name="html_file" label="${job_name}.html"/>
-
-</outputs>
-<configfiles>
-<configfile name="runMe">
-#homozygosity-mapping script
-#input file: Q200 variants no wildtype
-ourargs = commandArgs(TRUE)
-inf = ourargs[1]
-outf = ourargs[2]
-Q200var = read.table(inf)
-
-
-het &lt;- Q200var[grep("AF=0.5",Q200var\$V8), ]
-hom &lt;- Q200var[grep("AF=1.0",Q200var\$V8), ]
-
-get_hom1 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr1")
-  het &lt;- subset(x2,x2\$V1=="chr1")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom2 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr2")
-  het &lt;- subset(x2,x2\$V1=="chr2")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom3 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr3")
-  het &lt;- subset(x2,x2\$V1=="chr3")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom4 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr4")
-  het &lt;- subset(x2,x2\$V1=="chr4")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom5 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr5")
-  het &lt;- subset(x2,x2\$V1=="chr5")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom6 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr6")
-  het &lt;- subset(x2,x2\$V1=="chr6")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom7 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr7")
-  het &lt;- subset(x2,x2\$V1=="chr7")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom8 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr8")
-  het &lt;- subset(x2,x2\$V1=="chr8")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom9 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr9")
-  het &lt;- subset(x2,x2\$V1=="chr9")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom10 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr10")
-  het &lt;- subset(x2,x2\$V1=="chr10")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom11 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr11")
-  het &lt;- subset(x2,x2\$V1=="chr11")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom12 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr12")
-  het &lt;- subset(x2,x2\$V1=="chr12")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom13 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr13")
-  het &lt;- subset(x2,x2\$V1=="chr13")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom14 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr14")
-  het &lt;- subset(x2,x2\$V1=="chr14")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom15 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr15")
-  het &lt;- subset(x2,x2\$V1=="chr15")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom16 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr16")
-  het &lt;- subset(x2,x2\$V1=="chr16")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom17 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr17")
-  het &lt;- subset(x2,x2\$V1=="chr17")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom18 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr18")
-  het &lt;- subset(x2,x2\$V1=="chr18")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom19 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr19")
-  het &lt;- subset(x2,x2\$V1=="chr19")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom20 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr20")
-  het &lt;- subset(x2,x2\$V1=="chr20")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom21 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr21")
-  het &lt;- subset(x2,x2\$V1=="chr21")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom22 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr22")
-  het &lt;- subset(x2,x2\$V1=="chr22")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom23 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr23")
-  het &lt;- subset(x2,x2\$V1=="chr23")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom24 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr24")
-  het &lt;- subset(x2,x2\$V1=="chr24")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-get_hom25 &lt;- function(x1,x2){
-  hom &lt;- subset(x1,x1\$V1=="chr25")
-  het &lt;- subset(x2,x2\$V1=="chr25")
-  d1 &lt;- density(hom\$V2)
-  d2 &lt;- density(het\$V2)
-  return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
-}
-
-h1 &lt;- get_hom1(hom,het)
-h2 &lt;- get_hom2(hom,het)
-h3 &lt;- get_hom3(hom,het)
-h4 &lt;- get_hom4(hom,het)
-h5 &lt;- get_hom5(hom,het)
-h6 &lt;- get_hom6(hom,het)
-h7 &lt;- get_hom7(hom,het)
-h8 &lt;- get_hom8(hom,het)
-h9 &lt;- get_hom9(hom,het)
-h10 &lt;- get_hom10(hom,het)
-h11 &lt;- get_hom11(hom,het)
-h12 &lt;- get_hom12(hom,het)
-h13 &lt;- get_hom13(hom,het)
-h14 &lt;- get_hom14(hom,het)
-h15 &lt;- get_hom15(hom,het)
-h16 &lt;- get_hom16(hom,het)
-h17 &lt;- get_hom17(hom,het)
-h18 &lt;- get_hom18(hom,het)
-h19 &lt;- get_hom19(hom,het)
-h20 &lt;- get_hom20(hom,het)
-h21 &lt;- get_hom21(hom,het)
-h22 &lt;- get_hom22(hom,het)
-h23 &lt;- get_hom23(hom,het)
-h24 &lt;- get_hom24(hom,het)
-h25 &lt;- get_hom25(hom,het)
-
-
-
-#Extracting maxima (homozygosity-mapping)
-
-zero_crossings_maxima &lt;- function(x,y){
-
-  dy &lt;- diff(y)
-  xn &lt;- x[2:length(x)]
-
-  n &lt;- length(dy)
-  t1 &lt;- dy[1:n-1];
-  t2 &lt;- dy[2:n];
-  tt &lt;- t1*t2
-  ind &lt;- which(tt &lt; 0)
-
-  # print out the maxima
-  d2y &lt;- diff(dy)
-  vind &lt;- c()
-  vy &lt;- c()
-  vxn &lt;- c()
-
-  xnn &lt;- x[3:length(x)]
-  for(i in 1:length(ind)){
-    ### cat( ind[i], d2y[ind[i]], '\n')
-    if( d2y[ind[i]] &lt; 0 ){
-      ### cat(ind[i], y[ind[i]], xn[ind[i]], "\n" )
-      vind &lt;- c(vind, ind[i])
-      vy &lt;- c(vy, y[ind[i]])
-      vxn &lt;- c(vxn, xn[ind[i]])
-    }
-  }
-
-  maxima &lt;- data.frame(ind=vind, y=vy, x=vxn)
-  maxima &lt;- maxima[order(maxima\$y),]
-
-  #print( maxima )
-
-  results &lt;- list()
-  results\$dy &lt;- dy
-  results\$xn &lt;- xn
-  results\$maxima &lt;- maxima[ order(maxima\$y, decreasing=T), ]
-
-  return(results)
-}
-
-res1 &lt;- zero_crossings_maxima(h1\$x,h1\$ratio)
-res2 &lt;- zero_crossings_maxima(h2\$x,h2\$ratio)
-res3 &lt;- zero_crossings_maxima(h3\$x,h3\$ratio)
-res4 &lt;- zero_crossings_maxima(h4\$x,h4\$ratio)
-res5 &lt;- zero_crossings_maxima(h5\$x,h5\$ratio)
-res6 &lt;- zero_crossings_maxima(h6\$x,h6\$ratio)
-res7 &lt;- zero_crossings_maxima(h7\$x,h7\$ratio)
-res8 &lt;- zero_crossings_maxima(h8\$x,h8\$ratio)
-res9 &lt;- zero_crossings_maxima(h9\$x,h9\$ratio)
-res10 &lt;- zero_crossings_maxima(h10\$x,h10\$ratio)
-res11 &lt;- zero_crossings_maxima(h11\$x,h11\$ratio)
-res12 &lt;- zero_crossings_maxima(h12\$x,h12\$ratio)
-res13 &lt;- zero_crossings_maxima(h13\$x,h13\$ratio)
-res14 &lt;- zero_crossings_maxima(h14\$x,h14\$ratio)
-res15 &lt;- zero_crossings_maxima(h15\$x,h15\$ratio)
-res16 &lt;- zero_crossings_maxima(h16\$x,h16\$ratio)
-res17 &lt;- zero_crossings_maxima(h17\$x,h17\$ratio)
-res18 &lt;- zero_crossings_maxima(h18\$x,h18\$ratio)
-res19 &lt;- zero_crossings_maxima(h19\$x,h19\$ratio)
-res20 &lt;- zero_crossings_maxima(h20\$x,h20\$ratio)
-res21 &lt;- zero_crossings_maxima(h21\$x,h21\$ratio)
-res22 &lt;- zero_crossings_maxima(h22\$x,h22\$ratio)
-res23 &lt;- zero_crossings_maxima(h23\$x,h23\$ratio)
-res24 &lt;- zero_crossings_maxima(h24\$x,h24\$ratio)
-res25 &lt;- zero_crossings_maxima(h25\$x,h25\$ratio)
-
-
-#plotting (including maxima and adjusting to max ratio)
-# all on one plot
-pdf("allchrs.pdf",height=10, width=20)
-library(Hmisc)
-ymax &lt;- 1.05*max( c(h1\$ratio,h2\$ratio,h3\$ratio,h4\$ratio,h5\$ratio,h6\$ratio,h7\$ratio,h8\$ratio,h9\$ratio,h10\$ratio,h11\$ratio,h12\$ratio,h13\$ratio,h14\$ratio,h15\$ratio,h16\$ratio,h17\$ratio,h18\$ratio,h19\$ratio,h20\$ratio,h21\$ratio,h22\$ratio,h23\$ratio,h24\$ratio,h25\$ratio))
-par(mfrow=c(5,5))
-plot(h1\$x, h1\$ratio,type="l", lwd=3, col="red", xlim=c(0,61000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 1")
-points( res1\$maxima\$x[1], res1\$maxima\$y[1], col="red", pch=19 )
-abline(v=res1\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res1\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h2\$x, h2\$ratio,type="l", lwd=3, col="red", xlim=c(0,61000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 2")
-points( res2\$maxima\$x[1], res2\$maxima\$y[1], col="red", pch=19 )
-abline(v=res2\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res2\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h3\$x, h3\$ratio,type="l", lwd=3, col="red", xlim=c(0,64000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 3")
-points( res3\$maxima\$x[1], res3\$maxima\$y[1], col="red", pch=19 )
-abline(v=res3\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res3\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h4\$x, h4\$ratio,type="l", lwd=3, col="red",, xlim=c(0,64000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 4")
-points( res4\$maxima\$x[1], res4\$maxima\$y[1], col="red", pch=19 )
-abline(v=res4\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res4\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h5\$x, h5\$ratio,type="l", lwd=3, col="red",, xlim=c(0,76000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 5")
-points( res5\$maxima\$x[1], res5\$maxima\$y[1], col="red", pch=19 )
-abline(v=res5\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res5\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h6\$x, h6\$ratio,type="l", lwd=3, col="red",, xlim=c(0,60000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 6")
-points( res6\$maxima\$x[1], res6\$maxima\$y[1], col="red", pch=19 )
-abline(v=res6\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res6\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h7\$x, h7\$ratio,type="l", lwd=3, col="red",, xlim=c(0,78000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 7")
-points( res7\$maxima\$x[1], res7\$maxima\$y[1], col="red", pch=19 )
-abline(v=res7\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res7\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h8\$x, h8\$ratio,type="l", lwd=3, col="red",, xlim=c(0,57000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 8")
-points( res8\$maxima\$x[1], res8\$maxima\$y[1], col="red", pch=19 )
-abline(v=res8\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res8\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h9\$x, h9\$ratio,type="l", lwd=3, col="red",, xlim=c(0,59000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 9")
-points( res9\$maxima\$x[1], res9\$maxima\$y[1], col="red", pch=19 )
-abline(v=res9\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res9\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h10\$x, h10\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 10")
-points( res10\$maxima\$x[1], res10\$maxima\$y[1], col="red", pch=19 )
-abline(v=res10\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res10\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h11\$x, h11\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 11")
-points( res11\$maxima\$x[1], res11\$maxima\$y[1], col="red", pch=19 )
-abline(v=res11\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res11\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h12\$x, h12\$ratio,type="l", lwd=3, col="red",, xlim=c(0,51000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 12")
-points( res12\$maxima\$x[1], res12\$maxima\$y[1], col="red", pch=19 )
-abline(v=res12\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res12\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h13\$x, h13\$ratio,type="l", lwd=3, col="red",, xlim=c(0,55000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 13")
-points( res13\$maxima\$x[1], res13\$maxima\$y[1], col="red", pch=19 )
-abline(v=res13\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res13\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h14\$x, h14\$ratio,type="l", lwd=3, col="red",, xlim=c(0,54000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 14")
-points( res14\$maxima\$x[1], res14\$maxima\$y[1], col="red", pch=19 )
-abline(v=res14\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res14\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h15\$x, h15\$ratio,type="l", lwd=3, col="red",, xlim=c(0,48000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 15")
-points( res15\$maxima\$x[1], res15\$maxima\$y[1], col="red", pch=19 )
-abline(v=res15\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res15\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h16\$x, h16\$ratio,type="l", lwd=3, col="red",, xlim=c(0,59000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 16")
-points( res16\$maxima\$x[1], res16\$maxima\$y[1], col="red", pch=19 )
-abline(v=res16\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res16\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h17\$x, h17\$ratio,type="l", lwd=3, col="red",, xlim=c(0,54000000),, ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 17")
-points( res17\$maxima\$x[1], res17\$maxima\$y[1], col="red", pch=19 )
-abline(v=res17\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res17\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h18\$x, h18\$ratio,type="l", lwd=3, col="red", xlim=c(0,50000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 18")
-points( res18\$maxima\$x[1], res18\$maxima\$y[1], col="red", pch=19 )
-abline(v=res18\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res18\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h19\$x, h19\$ratio,type="l", lwd=3, col="red",, xlim=c(0,51000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", , main="Chromosome 19")
-points( res19\$maxima\$x[1], res19\$maxima\$y[1], col="red", pch=19 )
-abline(v=res19\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res19\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h20\$x, h20\$ratio,type="l", lwd=3, col="red",, xlim=c(0,56000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 20")
-points( res20\$maxima\$x[1], res20\$maxima\$y[1], col="red", pch=19 )
-abline(v=res20\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res20\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h21\$x, h21\$ratio,type="l", lwd=3, col="red",, xlim=c(0,45000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 21")
-points( res21\$maxima\$x[1], res21\$maxima\$y[1], col="red", pch=19 )
-abline(v=res21\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res21\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h22\$x, h22\$ratio,type="l", lwd=3, col="red",, xlim=c(0,43000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 22")
-points( res22\$maxima\$x[1], res22\$maxima\$y[1], col="red", pch=19 )
-abline(v=res22\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res22\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h23\$x, h23\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 23")
-points( res23\$maxima\$x[1], res23\$maxima\$y[1], col="red", pch=19 )
-abline(v=res23\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res23\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h24\$x, h24\$ratio,type="l", lwd=3, col="red",, xlim=c(0,44000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 24")
-points( res24\$maxima\$x[1], res24\$maxima\$y[1], col="red", pch=19 )
-abline(v=res24\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res24\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-plot(h25\$x, h25\$ratio,type="l", lwd=3, col="red",, xlim=c(0,39000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 25")
-points( res25\$maxima\$x[1], res25\$maxima\$y[1], col="red", pch=19 )
-abline(v=res25\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n",  paste(format(round(res25\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(10)
-dev.off()
-
-# seperate plots
-pdf("chr1.pdf",height=5, width=10)
-plot(h1\$x, h1\$ratio,type="l", lwd=3, col="red", xlim=c(0,61000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 1")
-points( res1\$maxima\$x[1], res1\$maxima\$y[1], col="red", pch=19 )
-abline(v=res1\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res1\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr2.pdf",height=5, width=10)
-plot(h2\$x, h2\$ratio,type="l", lwd=3, col="red", xlim=c(0,61000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 2")
-points( res2\$maxima\$x[1], res2\$maxima\$y[1], col="red", pch=19 )
-abline(v=res2\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res2\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr3.pdf",height=5, width=10)
-plot(h3\$x, h3\$ratio,type="l", lwd=3, col="red", xlim=c(0,64000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 3")
-points( res3\$maxima\$x[1], res3\$maxima\$y[1], col="red", pch=19 )
-abline(v=res3\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res3\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr4.pdf",height=5, width=10)
-plot(h4\$x, h4\$ratio,type="l", lwd=3, col="red",, xlim=c(0,64000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 4")
-points( res4\$maxima\$x[1], res4\$maxima\$y[1], col="red", pch=19 )
-abline(v=res4\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res4\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr5.pdf",height=5, width=10)
-plot(h5\$x, h5\$ratio,type="l", lwd=3, col="red",, xlim=c(0,76000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 5")
-points( res5\$maxima\$x[1], res5\$maxima\$y[1], col="red", pch=19 )
-abline(v=res5\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res5\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr6.pdf",height=5, width=10)
-plot(h6\$x, h6\$ratio,type="l", lwd=3, col="red",, xlim=c(0,60000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 6")
-points( res6\$maxima\$x[1], res6\$maxima\$y[1], col="red", pch=19 )
-abline(v=res6\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res6\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr7.pdf",height=5, width=10)
-plot(h7\$x, h7\$ratio,type="l", lwd=3, col="red",, xlim=c(0,78000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 7")
-points( res7\$maxima\$x[1], res7\$maxima\$y[1], col="red", pch=19 )
-abline(v=res7\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res7\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr8.pdf",height=5, width=10)
-plot(h8\$x, h8\$ratio,type="l", lwd=3, col="red",, xlim=c(0,57000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 8")
-points( res8\$maxima\$x[1], res8\$maxima\$y[1], col="red", pch=19 )
-abline(v=res8\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res8\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr9.pdf",height=5, width=10)
-plot(h9\$x, h9\$ratio,type="l", lwd=3, col="red",, xlim=c(0,59000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 9")
-points( res9\$maxima\$x[1], res9\$maxima\$y[1], col="red", pch=19 )
-abline(v=res9\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res9\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr10.pdf",height=5, width=10)
-plot(h10\$x, h10\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 10")
-points( res10\$maxima\$x[1], res10\$maxima\$y[1], col="red", pch=19 )
-abline(v=res10\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res10\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr11.pdf",height=5, width=10)
-plot(h11\$x, h11\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 11")
-points( res11\$maxima\$x[1], res11\$maxima\$y[1], col="red", pch=19 )
-abline(v=res11\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res11\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr12.pdf",height=5, width=10)
-plot(h12\$x, h12\$ratio,type="l", lwd=3, col="red",, xlim=c(0,51000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 12")
-points( res12\$maxima\$x[1], res12\$maxima\$y[1], col="red", pch=19 )
-abline(v=res12\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res12\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr13.pdf",height=5, width=10)
-plot(h13\$x, h13\$ratio,type="l", lwd=3, col="red",, xlim=c(0,55000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 13")
-points( res13\$maxima\$x[1], res13\$maxima\$y[1], col="red", pch=19 )
-abline(v=res13\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res13\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr14.pdf",height=5, width=10)
-plot(h14\$x, h14\$ratio,type="l", lwd=3, col="red",, xlim=c(0,54000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 14")
-points( res14\$maxima\$x[1], res14\$maxima\$y[1], col="red", pch=19 )
-abline(v=res14\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res14\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr15.pdf",height=5, width=10)
-plot(h15\$x, h15\$ratio,type="l", lwd=3, col="red",, xlim=c(0,48000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 15")
-points( res15\$maxima\$x[1], res15\$maxima\$y[1], col="red", pch=19 )
-abline(v=res15\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res15\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr16.pdf",height=5, width=10)
-plot(h16\$x, h16\$ratio,type="l", lwd=3, col="red",, xlim=c(0,59000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 16")
-points( res16\$maxima\$x[1], res16\$maxima\$y[1], col="red", pch=19 )
-abline(v=res16\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res16\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr17.pdf",height=5, width=10)
-plot(h17\$x, h17\$ratio,type="l", lwd=3, col="red",, xlim=c(0,54000000),, ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 17")
-points( res17\$maxima\$x[1], res17\$maxima\$y[1], col="red", pch=19 )
-abline(v=res17\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res17\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr18.pdf",height=5, width=10)
-plot(h18\$x, h18\$ratio,type="l", lwd=3, col="red", xlim=c(0,50000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 18")
-points( res18\$maxima\$x[1], res18\$maxima\$y[1], col="red", pch=19 )
-abline(v=res18\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res18\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr19.pdf",height=5, width=10)
-plot(h19\$x, h19\$ratio,type="l", lwd=3, col="red",, xlim=c(0,51000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", , main="Chromosome 19")
-points( res19\$maxima\$x[1], res19\$maxima\$y[1], col="red", pch=19 )
-abline(v=res19\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res19\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr20.pdf",height=5, width=10)
-plot(h20\$x, h20\$ratio,type="l", lwd=3, col="red",, xlim=c(0,56000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 20")
-points( res20\$maxima\$x[1], res20\$maxima\$y[1], col="red", pch=19 )
-abline(v=res20\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res20\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr21.pdf",height=5, width=10)
-plot(h21\$x, h21\$ratio,type="l", lwd=3, col="red",, xlim=c(0,45000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 21")
-points( res21\$maxima\$x[1], res21\$maxima\$y[1], col="red", pch=19 )
-abline(v=res21\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res21\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr22.pdf",height=5, width=10)
-plot(h22\$x, h22\$ratio,type="l", lwd=3, col="red",, xlim=c(0,43000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 22")
-points( res22\$maxima\$x[1], res22\$maxima\$y[1], col="red", pch=19 )
-abline(v=res22\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res22\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr23.pdf",height=5, width=10)
-plot(h23\$x, h23\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 23")
-points( res23\$maxima\$x[1], res23\$maxima\$y[1], col="red", pch=19 )
-abline(v=res23\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res23\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr24.pdf",height=5, width=10)
-plot(h24\$x, h24\$ratio,type="l", lwd=3, col="red",, xlim=c(0,44000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 24")
-points( res24\$maxima\$x[1], res24\$maxima\$y[1], col="red", pch=19 )
-abline(v=res24\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n", legend=paste(format(round(res24\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(nx=10)
-dev.off()
-pdf("chr25.pdf",height=5, width=10)
-plot(h25\$x, h25\$ratio,type="l", lwd=3, col="red",, xlim=c(0,39000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 25")
-points( res25\$maxima\$x[1], res25\$maxima\$y[1], col="red", pch=19 )
-abline(v=res25\$maxima\$x[1], lwd=1)
-legend("topleft", bty="n",  paste(format(round(res25\$maxima\$x[1] ,0),nsmall=0) ) )
-minor.tick(10)
-dev.off()
-</configfile>
-</configfiles>
-
-
-        <tests>
-        <test>
-        <param name="input1" value="cloudmap_zf_hom_mapping_test1_input.xls" ftype="vcf"/>
-        <param name="job_name" value="test1"/>
-        <param name="runMe" value="$runMe"/>
-        <output name="html_file" file="cloudmap_zf_hom_mapping_test1_output.html" ftype="html" lines_diff="5"/>
-        </test>
-        </tests>
-        
-
-<help>
-
-
-**What it Does**
-
-This tool is part of an update to the CloudMap pipeline (Minevich et al. 2012) to incorporate homozygosity mapping in zebrafish. This tool facilitates the mapping of zebrafish mutants from WGS or RNA-Seq data using bulk segregant homozygosity mapping. As input it takes a VCF list of variants (usually quality filtered for Q200 with WT SNPs/indels removed). This file is produced when running the CloudMap zebrafish pipeline. It plots an adjusted homozygosity ratio such that a high ratio indicates mapping region. It outputs an HTML file with combined and individual chromosome plots as PDFs. The adjusted homozygosity ratio is calculated as described in:
-
-IN SUBMISSION
-
-**Script**
-Pressing execute will run the following code over your input file and generate some outputs in your history::
-
-
- #homozygosity-mapping script
- #input file: Q200 variants no wildtype
- ourargs = commandArgs(TRUE)
- inf = ourargs[1]
- outf = ourargs[2]
- Q200var = read.table(inf)
- 
- 
- het &lt;- Q200var[grep("AF=0.5",Q200var\$V8), ]
- hom &lt;- Q200var[grep("AF=1.0",Q200var\$V8), ]
- 
- get_hom1 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr1")
-   het &lt;- subset(x2,x2\$V1=="chr1")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom2 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr2")
-   het &lt;- subset(x2,x2\$V1=="chr2")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom3 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr3")
-   het &lt;- subset(x2,x2\$V1=="chr3")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom4 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr4")
-   het &lt;- subset(x2,x2\$V1=="chr4")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom5 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr5")
-   het &lt;- subset(x2,x2\$V1=="chr5")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom6 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr6")
-   het &lt;- subset(x2,x2\$V1=="chr6")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom7 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr7")
-   het &lt;- subset(x2,x2\$V1=="chr7")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom8 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr8")
-   het &lt;- subset(x2,x2\$V1=="chr8")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom9 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr9")
-   het &lt;- subset(x2,x2\$V1=="chr9")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom10 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr10")
-   het &lt;- subset(x2,x2\$V1=="chr10")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom11 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr11")
-   het &lt;- subset(x2,x2\$V1=="chr11")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom12 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr12")
-   het &lt;- subset(x2,x2\$V1=="chr12")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom13 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr13")
-   het &lt;- subset(x2,x2\$V1=="chr13")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom14 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr14")
-   het &lt;- subset(x2,x2\$V1=="chr14")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom15 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr15")
-   het &lt;- subset(x2,x2\$V1=="chr15")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom16 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr16")
-   het &lt;- subset(x2,x2\$V1=="chr16")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom17 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr17")
-   het &lt;- subset(x2,x2\$V1=="chr17")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom18 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr18")
-   het &lt;- subset(x2,x2\$V1=="chr18")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom19 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr19")
-   het &lt;- subset(x2,x2\$V1=="chr19")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom20 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr20")
-   het &lt;- subset(x2,x2\$V1=="chr20")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom21 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr21")
-   het &lt;- subset(x2,x2\$V1=="chr21")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom22 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr22")
-   het &lt;- subset(x2,x2\$V1=="chr22")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom23 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr23")
-   het &lt;- subset(x2,x2\$V1=="chr23")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom24 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr24")
-   het &lt;- subset(x2,x2\$V1=="chr24")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- get_hom25 &lt;- function(x1,x2){
-   hom &lt;- subset(x1,x1\$V1=="chr25")
-   het &lt;- subset(x2,x2\$V1=="chr25")
-   d1 &lt;- density(hom\$V2)
-   d2 &lt;- density(het\$V2)
-   return(list( x=d1\$x, ratio=(d1\$y/d2\$y)*(nrow(hom)/nrow(het))))
- }
- 
- h1 &lt;- get_hom1(hom,het)
- h2 &lt;- get_hom2(hom,het)
- h3 &lt;- get_hom3(hom,het)
- h4 &lt;- get_hom4(hom,het)
- h5 &lt;- get_hom5(hom,het)
- h6 &lt;- get_hom6(hom,het)
- h7 &lt;- get_hom7(hom,het)
- h8 &lt;- get_hom8(hom,het)
- h9 &lt;- get_hom9(hom,het)
- h10 &lt;- get_hom10(hom,het)
- h11 &lt;- get_hom11(hom,het)
- h12 &lt;- get_hom12(hom,het)
- h13 &lt;- get_hom13(hom,het)
- h14 &lt;- get_hom14(hom,het)
- h15 &lt;- get_hom15(hom,het)
- h16 &lt;- get_hom16(hom,het)
- h17 &lt;- get_hom17(hom,het)
- h18 &lt;- get_hom18(hom,het)
- h19 &lt;- get_hom19(hom,het)
- h20 &lt;- get_hom20(hom,het)
- h21 &lt;- get_hom21(hom,het)
- h22 &lt;- get_hom22(hom,het)
- h23 &lt;- get_hom23(hom,het)
- h24 &lt;- get_hom24(hom,het)
- h25 &lt;- get_hom25(hom,het)
- 
- 
- 
- #Extracting maxima (homozygosity-mapping)
- 
- zero_crossings_maxima &lt;- function(x,y){
- 
-   dy &lt;- diff(y)
-   xn &lt;- x[2:length(x)]
- 
-   n &lt;- length(dy)
-   t1 &lt;- dy[1:n-1];
-   t2 &lt;- dy[2:n];
-   tt &lt;- t1*t2
-   ind &lt;- which(tt &lt; 0)
- 
-   # print out the maxima
-   d2y &lt;- diff(dy)
-   vind &lt;- c()
-   vy &lt;- c()
-   vxn &lt;- c()
- 
-   xnn &lt;- x[3:length(x)]
-   for(i in 1:length(ind)){
-     ### cat( ind[i], d2y[ind[i]], '\n')
-     if( d2y[ind[i]] &lt; 0 ){
-       ### cat(ind[i], y[ind[i]], xn[ind[i]], "\n" )
-       vind &lt;- c(vind, ind[i])
-       vy &lt;- c(vy, y[ind[i]])
-       vxn &lt;- c(vxn, xn[ind[i]])
-     }
-   }
- 
-   maxima &lt;- data.frame(ind=vind, y=vy, x=vxn)
-   maxima &lt;- maxima[order(maxima\$y),]
- 
-   #print( maxima )
- 
-   results &lt;- list()
-   results\$dy &lt;- dy
-   results\$xn &lt;- xn
-   results\$maxima &lt;- maxima[ order(maxima\$y, decreasing=T), ]
- 
-   return(results)
- }
- 
- res1 &lt;- zero_crossings_maxima(h1\$x,h1\$ratio)
- res2 &lt;- zero_crossings_maxima(h2\$x,h2\$ratio)
- res3 &lt;- zero_crossings_maxima(h3\$x,h3\$ratio)
- res4 &lt;- zero_crossings_maxima(h4\$x,h4\$ratio)
- res5 &lt;- zero_crossings_maxima(h5\$x,h5\$ratio)
- res6 &lt;- zero_crossings_maxima(h6\$x,h6\$ratio)
- res7 &lt;- zero_crossings_maxima(h7\$x,h7\$ratio)
- res8 &lt;- zero_crossings_maxima(h8\$x,h8\$ratio)
- res9 &lt;- zero_crossings_maxima(h9\$x,h9\$ratio)
- res10 &lt;- zero_crossings_maxima(h10\$x,h10\$ratio)
- res11 &lt;- zero_crossings_maxima(h11\$x,h11\$ratio)
- res12 &lt;- zero_crossings_maxima(h12\$x,h12\$ratio)
- res13 &lt;- zero_crossings_maxima(h13\$x,h13\$ratio)
- res14 &lt;- zero_crossings_maxima(h14\$x,h14\$ratio)
- res15 &lt;- zero_crossings_maxima(h15\$x,h15\$ratio)
- res16 &lt;- zero_crossings_maxima(h16\$x,h16\$ratio)
- res17 &lt;- zero_crossings_maxima(h17\$x,h17\$ratio)
- res18 &lt;- zero_crossings_maxima(h18\$x,h18\$ratio)
- res19 &lt;- zero_crossings_maxima(h19\$x,h19\$ratio)
- res20 &lt;- zero_crossings_maxima(h20\$x,h20\$ratio)
- res21 &lt;- zero_crossings_maxima(h21\$x,h21\$ratio)
- res22 &lt;- zero_crossings_maxima(h22\$x,h22\$ratio)
- res23 &lt;- zero_crossings_maxima(h23\$x,h23\$ratio)
- res24 &lt;- zero_crossings_maxima(h24\$x,h24\$ratio)
- res25 &lt;- zero_crossings_maxima(h25\$x,h25\$ratio)
- 
- 
- #plotting (including maxima and adjusting to max ratio)
- # all on one plot
- pdf("allchrs.pdf",height=10, width=20)
- library(Hmisc)
- ymax &lt;- 1.05*max( c(h1\$ratio,h2\$ratio,h3\$ratio,h4\$ratio,h5\$ratio,h6\$ratio,h7\$ratio,h8\$ratio,h9\$ratio,h10\$ratio,h11\$ratio,h12\$ratio,h13\$ratio,h14\$ratio,h15\$ratio,h16\$ratio,h17\$ratio,h18\$ratio,h19\$ratio,h20\$ratio,h21\$ratio,h22\$ratio,h23\$ratio,h24\$ratio,h25\$ratio))
- par(mfrow=c(5,5))
- plot(h1\$x, h1\$ratio,type="l", lwd=3, col="red", xlim=c(0,61000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 1")
- points( res1\$maxima\$x[1], res1\$maxima\$y[1], col="red", pch=19 )
- abline(v=res1\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res1\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h2\$x, h2\$ratio,type="l", lwd=3, col="red", xlim=c(0,61000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 2")
- points( res2\$maxima\$x[1], res2\$maxima\$y[1], col="red", pch=19 )
- abline(v=res2\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res2\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h3\$x, h3\$ratio,type="l", lwd=3, col="red", xlim=c(0,64000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 3")
- points( res3\$maxima\$x[1], res3\$maxima\$y[1], col="red", pch=19 )
- abline(v=res3\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res3\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h4\$x, h4\$ratio,type="l", lwd=3, col="red",, xlim=c(0,64000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 4")
- points( res4\$maxima\$x[1], res4\$maxima\$y[1], col="red", pch=19 )
- abline(v=res4\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res4\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h5\$x, h5\$ratio,type="l", lwd=3, col="red",, xlim=c(0,76000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 5")
- points( res5\$maxima\$x[1], res5\$maxima\$y[1], col="red", pch=19 )
- abline(v=res5\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res5\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h6\$x, h6\$ratio,type="l", lwd=3, col="red",, xlim=c(0,60000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 6")
- points( res6\$maxima\$x[1], res6\$maxima\$y[1], col="red", pch=19 )
- abline(v=res6\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res6\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h7\$x, h7\$ratio,type="l", lwd=3, col="red",, xlim=c(0,78000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 7")
- points( res7\$maxima\$x[1], res7\$maxima\$y[1], col="red", pch=19 )
- abline(v=res7\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res7\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h8\$x, h8\$ratio,type="l", lwd=3, col="red",, xlim=c(0,57000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 8")
- points( res8\$maxima\$x[1], res8\$maxima\$y[1], col="red", pch=19 )
- abline(v=res8\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res8\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h9\$x, h9\$ratio,type="l", lwd=3, col="red",, xlim=c(0,59000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 9")
- points( res9\$maxima\$x[1], res9\$maxima\$y[1], col="red", pch=19 )
- abline(v=res9\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res9\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h10\$x, h10\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 10")
- points( res10\$maxima\$x[1], res10\$maxima\$y[1], col="red", pch=19 )
- abline(v=res10\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res10\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h11\$x, h11\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 11")
- points( res11\$maxima\$x[1], res11\$maxima\$y[1], col="red", pch=19 )
- abline(v=res11\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res11\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h12\$x, h12\$ratio,type="l", lwd=3, col="red",, xlim=c(0,51000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 12")
- points( res12\$maxima\$x[1], res12\$maxima\$y[1], col="red", pch=19 )
- abline(v=res12\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res12\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h13\$x, h13\$ratio,type="l", lwd=3, col="red",, xlim=c(0,55000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 13")
- points( res13\$maxima\$x[1], res13\$maxima\$y[1], col="red", pch=19 )
- abline(v=res13\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res13\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h14\$x, h14\$ratio,type="l", lwd=3, col="red",, xlim=c(0,54000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 14")
- points( res14\$maxima\$x[1], res14\$maxima\$y[1], col="red", pch=19 )
- abline(v=res14\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res14\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h15\$x, h15\$ratio,type="l", lwd=3, col="red",, xlim=c(0,48000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 15")
- points( res15\$maxima\$x[1], res15\$maxima\$y[1], col="red", pch=19 )
- abline(v=res15\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res15\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h16\$x, h16\$ratio,type="l", lwd=3, col="red",, xlim=c(0,59000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 16")
- points( res16\$maxima\$x[1], res16\$maxima\$y[1], col="red", pch=19 )
- abline(v=res16\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res16\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h17\$x, h17\$ratio,type="l", lwd=3, col="red",, xlim=c(0,54000000),, ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 17")
- points( res17\$maxima\$x[1], res17\$maxima\$y[1], col="red", pch=19 )
- abline(v=res17\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res17\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h18\$x, h18\$ratio,type="l", lwd=3, col="red", xlim=c(0,50000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 18")
- points( res18\$maxima\$x[1], res18\$maxima\$y[1], col="red", pch=19 )
- abline(v=res18\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res18\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h19\$x, h19\$ratio,type="l", lwd=3, col="red",, xlim=c(0,51000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", , main="Chromosome 19")
- points( res19\$maxima\$x[1], res19\$maxima\$y[1], col="red", pch=19 )
- abline(v=res19\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res19\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h20\$x, h20\$ratio,type="l", lwd=3, col="red",, xlim=c(0,56000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 20")
- points( res20\$maxima\$x[1], res20\$maxima\$y[1], col="red", pch=19 )
- abline(v=res20\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res20\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h21\$x, h21\$ratio,type="l", lwd=3, col="red",, xlim=c(0,45000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 21")
- points( res21\$maxima\$x[1], res21\$maxima\$y[1], col="red", pch=19 )
- abline(v=res21\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res21\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h22\$x, h22\$ratio,type="l", lwd=3, col="red",, xlim=c(0,43000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 22")
- points( res22\$maxima\$x[1], res22\$maxima\$y[1], col="red", pch=19 )
- abline(v=res22\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res22\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h23\$x, h23\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 23")
- points( res23\$maxima\$x[1], res23\$maxima\$y[1], col="red", pch=19 )
- abline(v=res23\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res23\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h24\$x, h24\$ratio,type="l", lwd=3, col="red",, xlim=c(0,44000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 24")
- points( res24\$maxima\$x[1], res24\$maxima\$y[1], col="red", pch=19 )
- abline(v=res24\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res24\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- plot(h25\$x, h25\$ratio,type="l", lwd=3, col="red",, xlim=c(0,39000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 25")
- points( res25\$maxima\$x[1], res25\$maxima\$y[1], col="red", pch=19 )
- abline(v=res25\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n",  paste(format(round(res25\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(10)
- dev.off()
- 
- # seperate plots
- pdf("chr1.pdf",height=5, width=10)
- plot(h1\$x, h1\$ratio,type="l", lwd=3, col="red", xlim=c(0,61000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 1")
- points( res1\$maxima\$x[1], res1\$maxima\$y[1], col="red", pch=19 )
- abline(v=res1\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res1\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr2.pdf",height=5, width=10)
- plot(h2\$x, h2\$ratio,type="l", lwd=3, col="red", xlim=c(0,61000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 2")
- points( res2\$maxima\$x[1], res2\$maxima\$y[1], col="red", pch=19 )
- abline(v=res2\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res2\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr3.pdf",height=5, width=10)
- plot(h3\$x, h3\$ratio,type="l", lwd=3, col="red", xlim=c(0,64000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 3")
- points( res3\$maxima\$x[1], res3\$maxima\$y[1], col="red", pch=19 )
- abline(v=res3\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res3\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr4.pdf",height=5, width=10)
- plot(h4\$x, h4\$ratio,type="l", lwd=3, col="red",, xlim=c(0,64000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 4")
- points( res4\$maxima\$x[1], res4\$maxima\$y[1], col="red", pch=19 )
- abline(v=res4\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res4\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr5.pdf",height=5, width=10)
- plot(h5\$x, h5\$ratio,type="l", lwd=3, col="red",, xlim=c(0,76000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 5")
- points( res5\$maxima\$x[1], res5\$maxima\$y[1], col="red", pch=19 )
- abline(v=res5\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res5\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr6.pdf",height=5, width=10)
- plot(h6\$x, h6\$ratio,type="l", lwd=3, col="red",, xlim=c(0,60000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 6")
- points( res6\$maxima\$x[1], res6\$maxima\$y[1], col="red", pch=19 )
- abline(v=res6\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res6\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr7.pdf",height=5, width=10)
- plot(h7\$x, h7\$ratio,type="l", lwd=3, col="red",, xlim=c(0,78000000),ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 7")
- points( res7\$maxima\$x[1], res7\$maxima\$y[1], col="red", pch=19 )
- abline(v=res7\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res7\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr8.pdf",height=5, width=10)
- plot(h8\$x, h8\$ratio,type="l", lwd=3, col="red",, xlim=c(0,57000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 8")
- points( res8\$maxima\$x[1], res8\$maxima\$y[1], col="red", pch=19 )
- abline(v=res8\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res8\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr9.pdf",height=5, width=10)
- plot(h9\$x, h9\$ratio,type="l", lwd=3, col="red",, xlim=c(0,59000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio" , main="Chromosome 9")
- points( res9\$maxima\$x[1], res9\$maxima\$y[1], col="red", pch=19 )
- abline(v=res9\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res9\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr10.pdf",height=5, width=10)
- plot(h10\$x, h10\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 10")
- points( res10\$maxima\$x[1], res10\$maxima\$y[1], col="red", pch=19 )
- abline(v=res10\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res10\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr11.pdf",height=5, width=10)
- plot(h11\$x, h11\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 11")
- points( res11\$maxima\$x[1], res11\$maxima\$y[1], col="red", pch=19 )
- abline(v=res11\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res11\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr12.pdf",height=5, width=10)
- plot(h12\$x, h12\$ratio,type="l", lwd=3, col="red",, xlim=c(0,51000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 12")
- points( res12\$maxima\$x[1], res12\$maxima\$y[1], col="red", pch=19 )
- abline(v=res12\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res12\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr13.pdf",height=5, width=10)
- plot(h13\$x, h13\$ratio,type="l", lwd=3, col="red",, xlim=c(0,55000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 13")
- points( res13\$maxima\$x[1], res13\$maxima\$y[1], col="red", pch=19 )
- abline(v=res13\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res13\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr14.pdf",height=5, width=10)
- plot(h14\$x, h14\$ratio,type="l", lwd=3, col="red",, xlim=c(0,54000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 14")
- points( res14\$maxima\$x[1], res14\$maxima\$y[1], col="red", pch=19 )
- abline(v=res14\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res14\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr15.pdf",height=5, width=10)
- plot(h15\$x, h15\$ratio,type="l", lwd=3, col="red",, xlim=c(0,48000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 15")
- points( res15\$maxima\$x[1], res15\$maxima\$y[1], col="red", pch=19 )
- abline(v=res15\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res15\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr16.pdf",height=5, width=10)
- plot(h16\$x, h16\$ratio,type="l", lwd=3, col="red",, xlim=c(0,59000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 16")
- points( res16\$maxima\$x[1], res16\$maxima\$y[1], col="red", pch=19 )
- abline(v=res16\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res16\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr17.pdf",height=5, width=10)
- plot(h17\$x, h17\$ratio,type="l", lwd=3, col="red",, xlim=c(0,54000000),, ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 17")
- points( res17\$maxima\$x[1], res17\$maxima\$y[1], col="red", pch=19 )
- abline(v=res17\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res17\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr18.pdf",height=5, width=10)
- plot(h18\$x, h18\$ratio,type="l", lwd=3, col="red", xlim=c(0,50000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 18")
- points( res18\$maxima\$x[1], res18\$maxima\$y[1], col="red", pch=19 )
- abline(v=res18\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res18\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr19.pdf",height=5, width=10)
- plot(h19\$x, h19\$ratio,type="l", lwd=3, col="red",, xlim=c(0,51000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", , main="Chromosome 19")
- points( res19\$maxima\$x[1], res19\$maxima\$y[1], col="red", pch=19 )
- abline(v=res19\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res19\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr20.pdf",height=5, width=10)
- plot(h20\$x, h20\$ratio,type="l", lwd=3, col="red",, xlim=c(0,56000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 20")
- points( res20\$maxima\$x[1], res20\$maxima\$y[1], col="red", pch=19 )
- abline(v=res20\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res20\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr21.pdf",height=5, width=10)
- plot(h21\$x, h21\$ratio,type="l", lwd=3, col="red",, xlim=c(0,45000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 21")
- points( res21\$maxima\$x[1], res21\$maxima\$y[1], col="red", pch=19 )
- abline(v=res21\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res21\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr22.pdf",height=5, width=10)
- plot(h22\$x, h22\$ratio,type="l", lwd=3, col="red",, xlim=c(0,43000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 22")
- points( res22\$maxima\$x[1], res22\$maxima\$y[1], col="red", pch=19 )
- abline(v=res22\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res22\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr23.pdf",height=5, width=10)
- plot(h23\$x, h23\$ratio,type="l", lwd=3, col="red",, xlim=c(0,47000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 23")
- points( res23\$maxima\$x[1], res23\$maxima\$y[1], col="red", pch=19 )
- abline(v=res23\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res23\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr24.pdf",height=5, width=10)
- plot(h24\$x, h24\$ratio,type="l", lwd=3, col="red",, xlim=c(0,44000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 24")
- points( res24\$maxima\$x[1], res24\$maxima\$y[1], col="red", pch=19 )
- abline(v=res24\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n", legend=paste(format(round(res24\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(nx=10)
- dev.off()
- pdf("chr25.pdf",height=5, width=10)
- plot(h25\$x, h25\$ratio,type="l", lwd=3, col="red",, xlim=c(0,39000000), ylim=c(0,ymax), xlab="Location (Mb)", ylab="Adjusted homozygosity ratio", main="Chromosome 25")
- points( res25\$maxima\$x[1], res25\$maxima\$y[1], col="red", pch=19 )
- abline(v=res25\$maxima\$x[1], lwd=1)
- legend("topleft", bty="n",  paste(format(round(res25\$maxima\$x[1] ,0),nsmall=0) ) )
- minor.tick(10)
- dev.off()
-
-**Attribution**
-This Galaxy tool was created by r.poole@ucl.ac.uk at 13/10/2014 13:41:11
-using the Galaxy Tool Factory.
-
-See https://bitbucket.org/fubar/galaxytoolfactory for details of that project
-Please cite: Creating re-usable tools from scripts: The Galaxy Tool Factory. Ross Lazarus; Antony Kaspi; Mark Ziemann; The Galaxy Team. 
-Bioinformatics 2012; doi: 10.1093/bioinformatics/bts573
-
-
-</help>
-<citations>
-    <citation type="doi">10.1534/genetics.112.144204</citation>
-    <citation type="doi">10.1093/bioinformatics/bts573</citation>
-</citations>
-</tool>