lh3 / psmc

Implementation of the Pairwise Sequentially Markovian Coalescent (PSMC) model
Other
146 stars 60 forks source link

Not issues and a question #30

Open dongyiyi opened 3 years ago

dongyiyi commented 3 years ago

Why does not the x-axis of PSMC plotting be divided into equal parts/length on the time scale line?

fredjaya commented 3 years ago

It's log scaled.

kishor2019 commented 3 years ago

Hi @fredjaya Would you like to help me to set up the generation time and mutation rate? ./PSMC.sh genome.fa num_of_thread R1.fq R2.fq run_folder sequence_depth psmc_PATH #this was my command line

This was my bash script

mkdir ${5}

1.create run folder, build index and alignment

bwa mem -t ${2} -o ${5}/psmcOut.sam ${1} ${3} ${4} samtools view -b -@ ${2} -o ${5}/psmcOut.bam ${5}/psmcOut.sam samtools sort -@ ${2} -o ${5}/psmcOut_sorted.bam ${5}/psmcOut.bam echo "1.create run folder, build index and alignment processed"

2.use samtools get diploid.fq.gz

samtools mpileup -C50 -uf ${1} ${5}/psmcOut_sorted.bam > ${5}/psmcOut.bcf bcftools call -c ${5}/psmcOut.bcf > ${5}/psmcOut.vcf depth2x=expr ${6} \* 2 vcfutils.pl vcf2fq -d ${7} -D ${depth2x} ${5}/psmcOut.vcf | gzip > ${5}/diploid.fq.gz echo "2.use samtools get diploid.fq.gz processed"

3.psmc

${7}/utils/fq2psmcfa -q20 ${5}/diploid.fq.gz > ${5}/diploid.psmcfa ${7}/psmc -N25 -t15 -r5 -p "4+25*2+4+6" -o ${5}/diploid.psmc ${5}/diploid.psmcfa ${7}/utils/psmc2history.pl ${5}/diploid.psmc | ${7}/utils/history2ms.pl > ${5}/ms-cmd.sh ${7}/utils/psmc_plot.pl ${5}/diploid ${5}/diploid.psmc echo "3.psmc processed!"

But i had received a result with g = 25, µ = 2.5*10E-8 and the image is not interpretable that differes from actual scenerio. As i am running for a fish, i want to change the generation time g =2 and µ = 2.22 × 10-9. How i can fix this? Thanks a lot. I hope i will receive a warm response.