Merge pull request #1 from pdimens/anchor

Almost complete rewrite

.gitignore

@@ -1,6 +1,20 @@
 /.snakemake
-/maps.splitchrom
+.vscode
+/1_ParentCall
+/2_Filtering
+/3_SeparateChromosomes
+/4_OrderMarkers
+/5_Trim
+/6_OrderMarkers
+/7_Distances
+/7_DistancesSexAverage
+/7_Intervals
+/8_RepeatMask
+/9_Chain
+/10_Anchoring
+/11_MareyMaps
+/JoinSingles2All_iter
+/RefineMap
 /pedigree.txt
-/data.call.gz
-/data_f.call.gz
+map.master
 /*.vcf

LepWrap

@@ -10,9 +10,18 @@ fi
 if [[ -z "$1" ]]; then
   echo "Perform the modules of Lep-Map3. Make sure config.yaml is properly configured for how you intend to run things."
   echo ""
-  echo "[usage] LepMWrap <number of threads to use>"
+  echo "[usage] LepWrap <number of threads to use>"
   echo "[example] LepWrap 16" 
   exit 1
 fi
 
-snakemake --cores $1 --snakefile ./rules/Snakefile --directory .
+echo "Running Lep-Map3"
+sleep 2s
+snakemake --cores $1 --snakefile ./rules/LepMap3.smk --directory .
+
+LA=$(grep "run_lepanchor" config.yaml | cut -d":" -f2 | xargs | tr '[:upper:]' '[:lower:]')
+if [ $LA == "true" ]; then
+echo -e "\nRunning Lep-Anchor"
+sleep 2s
+snakemake --cores $1 --snakefile ./rules/LepAnchor.smk --directory .
+fi

config.yaml

@@ -1,20 +1,24 @@
 # Configuration file for LepWrap
 
+###################
+#### Lep-Map 3 ####
+###################
+
 #---------------#
 #  ParentCall2  #
 #---------------#
 # the filtered VCF file with your genotype likelihoods:
-vcf: "out.14.missrecode.vcf"
+vcf: "tons_of_snps.vcf"
 
-# Instructions to create pedigree file: https://sourceforge.net/p/lep-map3/wiki/LM3 Home/#parentcall2
+# Instructions to create pedigree file: https://sourceforge.net/p/lep-map3/wiki/software/LepMap3 Home/#parentcall2
 # the pedigree file associated with your data
 pedigree: "pedigree.txt"
 
 #---------------#
 #   Filtering2  #
 #---------------#
-# set this to 0 if you want to skip the Filtering2 module
-data_tol: 0.01
+# data tolerance value-- set this to 0 if you want to skip the Filtering2 module
+data_tol: 0.1
 
 #------------------------#
 #  SeperateChromosomes2  #
@@ -23,17 +27,20 @@ data_tol: 0.01
 # LOD score in the range of lod_min to lod_max
 
 # The minimum LOD for SeperateChromosomes2
-lod_min: 1
+lod_min: 5
 
 # The maximum LOD for SeperateChromosomes2
-lod_max: 40
+lod_max: 30
 
 # Use only markers with informative father (1), mother(2), both parents(3) or neither parent(0)
 informative: "informativeMask=123"
 
 #-------------------#
 #  JoinSingles2ALL  #
 #-------------------#
+# set this to false if you want to skip joining singles (0) to linkage groups
+run_joinsingles2all: true
+
 # these are the parameters for JoinSingles2ALL, and are highly data-dependent
 # start with lower values for lod_limit and increase as necessary
 lod_limit: "lodLimit=2"
@@ -49,16 +56,20 @@ lod_difference: "lodDifference=2"
 # Set exp_lg to your expected number of chromosomes
 exp_lg: 24
 
-# Set iter to the number of iterations you want per chromosome (more is better). Recommend 30 or more
-iter: 100
+# Set iterations to the number of iterations you want per chromosome (more is better). Recommend 30 or more
+iterations: 5
 
-# Set threads_per to the number of threads you would like to use per iteration
-threads_per: 2
+# Set threads_per to the number of threads you would like to use per linkage group for ordering
+threads_per: 5
 
 # Choose your distance method by commenting the line you dont want
 dist_method: "useKosambi=1"
 #dist_method: "useMorgan=1"
 
+# number of iterations to use of OrderMarkers2 phasing.
+# this number is typically 1-2
+phase_iterations: 1
+
 #-----------------#
 #  Edge Trimming  #
 #-----------------#
@@ -71,3 +82,33 @@ edge_length: 15
 
 # Set trim_cuttoff to the centiMorgan distance cutoff (10 is reasonable)
 trim_cutoff: 10 
+
+
+
+####################
+#### Lep-Anchor ####
+####################
+# change this to true if you also want to run Lep-Anchor
+run_lepanchor: false
+
+# the path to the genome assembly you are trying to anchor
+# ideally it *is not* gzipped and ends in .fa, but there are built-in workarounds if it is. 
+assembly: "boa.constrictor.fasta"
+
+# the number of linkage groups you have
+lg_count: 24
+
+# if you have a PAF file of long reads mapped to your genome, add it here, otherwise leave the text as "/dev/null"
+PAF_file: "alnPB.paf"
+
+# if you have a proximity file add it here, otherwise leave the text as "/dev/null".
+# This isn't yet implemented in the official Lep-Anchor wrapper, so don't change this.
+proximity_file: "/dev/null"
+
+# Which system is appropriate for running the HaploMerger2 section? Comment out the one that doesn't apply.
+# This is kinda necessary for the correct lastz binaries. If you have a Debian or Arch based system, Ubuntu
+# should work (tested for Arch Linux, at least). If using a Fedora/RedHat based system, try the CentOS options.
+# If you aren't sure, try "Ubuntu" and see what happens! I'm not actually convinced it matters, but you never know.
+OS_info: "Ubuntu"
+#OS_info: "CentOS5"
+#OS_info: "CentOS6"

rules/LAtest.smk

@@ -0,0 +1,9 @@
+rule validation:
+  input:
+  output:
+  message:
+  threads: 30
+  shell:
+    """
+    awk -f software/LepAnchor/scripts/liftover.awk chr1.agp order1.input | sort -V | grep CHR > order1.liftover
+    """