Difference between revisions of "Tips kang"

Revision as of 07:42, 1 May 2014

Maker

$ /data2/k821209/programs/maker/bin/gff3_merge -d redbean.ref_master_datastore_index.log
$ cat redbean.ref.all.gff | awk '$2=="maker"' - > redbean.ref.all.gff.maker.gff3
$ python3 sort.py redbean.ref.all.gff.maker.gff3
$ python3 genename_change_nonAnchor.py
$ /data2/k821209/programs/maker/bin/map_gff_ids Vang.scaffold.map redbean.ref.all.gff

Orthomcl

orthomcl 는 알아서 이름을 바꾸면된다.
Make config file : orthomcl.config.template
dbVendor=mysql
dbConnectString=dbi:mysql:orthomcl:mysql_local_infile=1:localhost:3306
dbLogin=orthomcl
dbPassword=8804555
similarSequencesTable=SimilarSequences
orthologTable=Ortholog
inParalogTable=InParalog
coOrthologTable=CoOrtholog
interTaxonMatchView=InterTaxonMatch
percentMatchCutoff=50
evalueExponentCutoff=-5
oracleIndexTblSpc=NONE

$ mysql -u root -p
mysql> CREATE DATABASE orthomcl;
mysql> GRANT SELECT,INSERT,UPDATE,DELETE,CREATE VIEW,CREATE, INDEX, DROP on orthomcl.* TO orthomcl@localhost;
mysql> set password for orthomcl@localhost = password('8804555');
$ orthomclInstallSchema orthomcl.config.template
$ orthomclAdjustFasta VRA final.assembly.pep.fasta.primary.SH.fa 1
$ mkdir ./compliantFasta/
$ mv VRA.fasta ./compliantFasta/ # same for the ATH.fasta etc.
$ orthomclFilterFasta ./compliantFasta/ 10 20
$ blastall -m 8 -a 10 -p blastp -F 'm S' -v 100000 -b 100000 -z 259701 -e 1e-5 -d goodProteins.fasta -i goodProteins.fasta -o goodProteins.fasta.self.out # 27개 할 때 보름 이상걸림
$ orthomclBlastParser goodProteins.fasta.self.out compliantFasta/ >> similarSequences.txt # 27개 accession할 때 이 스텝만 6시간 정도 걸림
$ orthomclLoadBlast orthomcl.config.template similarSequences.txt
$ orthomclPairs orthomcl.config.template orthomcl_pairs.log cleanup=no
$ orthomclDumpPairsFiles orthomcl.config.template
$ mcl mclInput --abc -I 1.5 -o mclOutput
$ orthomclMclToGroups CLST 1000 < mclOutput > groups.txt

Errors
execute failed: Can't create/write to file /tmp 어쩌고..
 이것은 mysql 버전에러라고 하는군. my.cnf 에서 tmp 디렉토리를 data 디렉토리랑 이름을 같게 만들어주면 돌아감. 
Duplicate entry error
 $ sort -k3,4 similarSequences.txt.primary > similarSequences.txt.primary.sort
 $ python3 duplicates_remove_ver1.py # 63:/data2/k821209/good_python_code/ 
 현재 테스트중

Protest

193server$ cd /data2/k821209/programs/prottest-3.4-20140123/  # 이 디렉토리로 이동해서 실행. 라이브러리 위치 때문에 따로 지정하기 귀찮기 때문에.. 
193server$ parallel -j 10 java -jar prottest-3.4.jar -i {1} -o {1}.protest.out -all-distributions -F -AIC -BIC -tc 0.5  :::  /data2/k821209/Mungbean/Mungbean_assemble/ver6/mcscanx/beast_prepare/all/others_ver2/manual/*.nex

bamUtil

bam2fastq

bam bam2FastQ --in filename.bam

Picard

sam2fastq

java -jar SamToFastq.jar INPUT=[file_bam] INCLUDE_NON_PRIMARY_ALIGNMENTS=True FASTQ=Name_1.txt SECOND_END_FASTQ=Name_2.txt VALIDATION_STRINGENCY=LENIENT

BEAST

beauti

/data2/k821209/programs/BEAST/bin/beauti -template /data2/k821209/programs/BEAST/templates/StarBeast.xml

beast

'single locus'

java -Xms256m -Xmx202400m -Djava.library.path=/data/program/BEAST/lib/ -jar /data/program/BEAST/lib/beast.jar -beagle_SSE -beagle_instances 10 ver7.xml

'multi loci'

java -Xms256m -Xmx202400m -Djava.library.path=/data/program/BEAST/lib/ -jar /data/program/BEAST/lib/beast.jar -beagle_SSE ver7.xml
java -Xms256m -Xmx202400m -Djava.library.path=/data/program/BEAST/lib/ -jar /data/program/BEAST/lib/beast.jar -beagle_SSE 557.xml # 557 loci : Total calculation time: 373863.472 seconds 

GPU with multi loci

C:\Users\jun\Desktop\kangyangjae\BEAST.v2.1.2\BEAST\557>java -Xmx15g -jar ..\lib\beast.jar -beagle_order 0,0,0,0,0,0,0,1

0,0,0,0,0,0,0,1 의 의미는 "7개 loci정보 cpu에 던지고 8번째는 gpu에 던져라" loci가 많을 경우는 반복으로 할당. gpu에 다 던지면 좋겠지만 gpu메모리 한계 때문에 이렇게 해야함.

18개 loci, 27 species alignment 에 GPU memory 770 Mb 정도 소요.

NCBI에 WGS 올리기 팁

꼼수이다. 참고로. 일단 그냥 Whole genome seqeunce fasta file을 그냥 WGS submission portal로 넣어버리고 하루를 기다린다. 뭐가 잘못되었다고 구구절절 온다. 이렇게 하라고도 구구절절온다. 그렇게 해서 보내준다. 끝.

Mutation rate

Legumes

"If the older duplication is assumed to have occurred around 58Myr ago, then the calculated rate of silent mutations extending back to the duplication would be 5.17*10-3, similar to previous estimates of 5.2*10–3" ^[1]

Python

Fisher's exact test

from scipy import stats

oddsratio, pvalue = stats.fisher_exact([[A,B], [C, D]]) ^[2]

63:/home/k821209/py/NGS/vcfq2fa.py : vcfutil로 만들어진 fq 파일을 fa로 변환

Excel

=TEXT(2.2323,"(0.00)")

(2.23)

Softwares

SNP 반영된 fasta 만들기

1. java -jar /data/program/picard-tools-1.91/CreateSequenceDictionary.jar R=Gmax_189.fa O=Gmax_189.dict
   • 인덱싱 작업인듯
2. java -Xmx2g -jar /data/program/GenomeAnalysisTK-3.1-1/GenomeAnalysisTK.jar -R /data/ref/Gmax_189.fa -T FastaAlternateReferenceMaker -o JM14.bam.sort.bam.bcg.bcf.d2D50.vcf.fa --variant JM14.bam.sort.bam.bcg.bcf.d2D50.vcf
   • bcftools view [filename.bcf] | /data1/KimSue/Gmax_189/ver1/bcf/bcf2fa.py prefix # GATK의 FastaAlternateReferenceMaker 는 calling이 되지 않는 지역을 reference 서열로 가져오는 문제가 있음.

GATK pipe http://www.broadinstitute.org/gatk//events/2038/

1. bwa mem -M -t 10 Va.ref.fa ysp-2_1.fastq.gz ysp-2_2.fastq.gz | /data/program/samtools-0.1.19/samtools view -Sb - | /data/program/samtools-0.1.19/samtools sort - ysp.bwamem.Va.ref.fa.sort # GATK pipe는 -M 옵션이 필요
2. /data/program/jdk1.7.0_25/bin/java -jar /data/program/picard-tools-1.91/MarkDuplicates.jar INPUT=ysp.bwamem.Va.ref.fa.sort.bam OUTPUT=ysp.bwamem.Va.ref.fa.sort.bam.dedup.bam METRICS_FILE=metrics.txt MAX_FILE_HANDLES_FOR_READ_ENDS_MAP=1000 TMP_DIR=./tmp/
   • MarkDuplicates 는 PCR duplicates를 제거하기 위한 작업. 왜 제거해야하는지는 자세히는 모르지만 PCR bias 때문에 Variant calling statistical model 에 문제가 된다고 함. 특히 PCR중에 polymerase가 실수 했을 경우 실수한 자리가 PCR bias로 대량 생산되었을 경우 마치 SNP처럼 calling되버림. 제거해야한다는 것이 결론. 그러나 특별히 제작된 라이브러리에는 적용해서는 안된다. 예를들면 특정 사이트가 잘리게 만들어놓은 라이브러리? 류들은 PCR bias로 오해받을 가능성? 이 있는듯..
   • MAX_FILE_HANDLES_FOR_READ_ENDS_MAP=1000 값을 높여도 되는데 ulimit -n 값을 올려야함. 재부팅 해야된다해서 그냥 낮춰 쓰는중. 그렇게 느린줄 모르겠음. 높여서 안써봐서 그런가..
   • TMP_DIR=./tmp/ 이거 빼먹으면 용량 없음 크리를 먹는다.

Maker

1. /data2/k821209/programs/maker/bin/gff3_merge -d Va.ref_master_datastore_index.log
2. /data2/k821209/programs/maker/bin/maker_map_ids --prefix=Vang --iterate=1 --suffix=. Va.ref.all.gff > id_map.txt # 복잡하게 나오는 maker의 유전자 이름들을 심플하게 바꾸는 툴
   • cat final.assembly_40k_GSFLX_pseudo5k.fasta.nonATGC.all.gff | awk '$2=="maker"' - > final.assembly_40k_GSFLX_pseudo5k.fasta.nonATGC.all.gff.maker.gff
   • python3 gff_parse.py final.assembly_40k_GSFLX_pseudo5k.fasta.nonATGC.all.gff.maker.gff
   • python3 /data2/k821209/Redbean/maker_pseudo/Va.ref.maker.output/genename_change_nonAnchor.py # 이름이 맘에 안들게 바뀌어서 개인적으로 만든 툴
3. /data2/k821209/programs/maker/bin/map_gff_ids id_map.txt Va.ref.all.gff # 그리 나온 이름들을 gff 반영하는 툴
4. python3 /data2/k821209/Redbean/maker_pseudo/Va.ref.maker.output/header_change.py Va.ref.all.maker.proteins.fasta Vang.scaffold.map # 그리 나온 이름들을 fasta에 반영하는 툴
5. /data2/k821209/programs/maker/bin/iprscan2gff3 Va.ref.all.maker.proteins.fasta.tsv.hc.tsv Va.ref.all.gff > Va.ref.all.gff.ipr.gff # interpro result를 jbrowser에 들어가는 모양으로 만들어주는 툴

Deconseq ^[3]

1. Illumina read의 contamination을 확인한다.
2. 63:/data/program/deconseq-standalone-0.4.3
3. /usr/bin/perl deconseq.pl -keep_tmp_files -f 800_both.fq -dbs bact,vir,arch -dbs_retain gmax

ePCR

Re-PCR

1. $ famap -tN -b genome.famap ./superscaf.ver1.fa
2. $ fahash -b genome.hash -w 12 -f3 ${PWD}/genome.famap
3. Work> /data/program/e-PCR-2.3.12/re-PCR -S genome.hash -n1 -g1 SSR.sts -o SSR.sts.mapped

SSR.sts

Mungbean_SSR_ID_1 CAAAAACATGAGTTGCACACAA TCATAACGCAGAACAGCGAA

Mungbean_SSR_ID_2 ATGTGTGTGAGCACCTCGAC TTTGGCCATGCAAGATGTAA

Mungbean_SSR_ID_4 GCGGTTCACCTAGCCATAAA GGACCCTTCTGTGCGTGTAT

Mungbean_SSR_ID_5 GTTTGTGCTGCGGATTCTTT TTGGCAATTTGGACTAAGGC

Mungbean_SSR_ID_7 TTGACCCAAAACTTACCAATTT GCTAAGGACTGGGGGTCTTC

Mummer, alignment draft genome to finished genome

$nucmer --prefix=ref_qry ref.fasta qry.fasta

$show-coords -rcl ref_qry.delta > ref_qry.coords

$show-aligns ref_qry.delta refname qryname > ref_qry.aligns

$show-tiling ref_qry.delta > ref_qry.tiling

추천 논문

Paterson, A.H., M. Freeling, H. Tang and X. Wang. 2010. Insights from the comparison of plant genome sequences. Annual Review of Plant Biology 61: 349-372. ^[4]

1. ↑ doi:10.1038/nature08670
2. ↑ scipy, fisher's exact
3. ↑ deconseq
4. ↑ Paterson, A.H., M. Freeling, H. Tang and X. Wang. 2010. Insights from the comparison of plant genome sequences. Annual Review of Plant Biology 61: 349-372.