2021 8/20 help追加
MELTは、iiluminaのペアエンドデータを使いリファレンスに存在しないmobile elementを検出するツール。以前1000 genomeで使われていたが、その後バージョンアップにより様々なゲノムに対応するようになった。SGEの分散コンピュータ環境から、SGEを使わない環境でも動く高いスケーラビリティを備えている。依存するソフトウエアも最小限になっており、エラーが起こりにくい設計を特徴とするらしい。
原理としては、一般的な構造変化の検出法と同じものを用いているが、途中からトランスポゾン検出に特化したパイプラインになっている。流れは大きく４つに分けることができる。１、ペアリードのマッピング（--> <--）が異常なペア（discordant read pair）の探索。２、おかしペアリードをユーザーが与えたトランスポゾン配列にアライメント。３、クロモソーム上を"waiking"して、トランスポゾンの挿入位置を正確に求める。４、情報を集計し、トランスポゾンごとに出力。この２のアライメントでbowtie2を使うため、MELTをランするにはbowtie2が前もってインストールされている必要がある。
依存
ダウンロード
公式サイトのdownloadタブからダウンロードできる。
MELT - Home
解凍して中に入る。MELT.jarが実行ファイル。
> java -jar MELT.jar
MELTv2.2.2 - Perform transposon analysis.
java -jar MELT.jar <Runtime>
Possible options for <Runtime>:
BuildTransposonZIP Create an MEI.zip file for input into other MELT Runtimes.
Preprocess Process .bam files for MELT analysis.
Single Run a single .bam file through MELT analysis.
SGE Run multiple .bam files through MELT analysis using SGE.
IndivAnalysis Step 1 of MELT-SPLIT -- Discover Putative MEIs.
GroupAnalysis Step 2 of MELT-SPLIT -- Determine MEI Breakpoint and Information.
Genotype Step 3 of MELT-SPLIT -- Genotype Putative MEIs.
MakeVCF Step 4 of MELT-SPLIT -- Perform final filtering and output VCF.
Source Build a source L1 file for input to Transduction.
TransductionFind Preprocess bam files for L1 3' Transduction discovery in already run MELT data.
TransductionMerge Search for L1 3' Transductions in already run MELT data using preprocessed .trans files.
Deletion-Genotype Genotype a .bam file for reference MEI deletions.
Deletion-Merge Merge multiple reference deletion genotypes into final VCF.
CALU CAlu Alu analysis. Classify Alu according to family and subfamily.
LINEU LineU L1 analysis. Classify L1 according to family and subfamily.
GENERIC Determine differences between a consensus ME sequence and additional copies.
--help/-help/-h will print this message and exit
> java -jar /MELT.jar Single
$ java -jar /MELT.jar Single
Command Line:
MELT.jar Single
Start time: 2021/08/20 11:00:51
Performing MELT analysis...
Missing required options: bamfile, h, w, t, n
usage: java -jar MELT.jar Single <options>
MELTv2.2.2 - MELT-Single - Perform transposon analysis on a single sample.
-a Reads have been aligned with bwa-aln. [false]
-ac Remove ac0 sites from final VCF file. [false]
-b <arg> Exclusion list for chromosomes. A '/' seperated list: i.e. to exclude chromosomes 1,2, and 4, put -b 1/2/4. [null]
-bamfile <arg> Bam file for MEI analysis.
-bowtie <arg> Path to the bowtie2 binary if not in PATH [null].
-c <arg> Coverage level of supplied bam file. Overides internal MELT coverage calculation if set. [null]
-cov <arg> Standard deviation cutoff when calling final sites in integer format. [35]
-d <arg> Minumum length of chromosome/contig size for calling elements. [1000000]
-e <arg> Expected insert size between reads. [500]
-exome Run MELT in exome mode (raises initial sensitivity of IndivAnalysis step). [false]
-h <arg> Path to the reference sequence used to align reads.
-j <arg> Total percentage of sites allowed to be no call (in integer form i.e. 25 percent would be -i 25, not .25). [25]
-k BAM file(s) have already been processed for discordant pairs (suffixes .fq, .disc, and .disc.bai are already present for the bam file in -l). [false]
-mcmq <arg> Allow MELT to use MC/MQ tags at greater than X percentage prevalance in the original bam file. [95.0]
-n <arg> Path to the genome annotation.
-nocleanup Do not cleanup MELT intermediate files after running. [false]
-q Alignments are phred+64 quality encoding. [false]
-r <arg> Read length of the supplied bam file(s). [100]
-s <arg> Standard deviation cutoff for excluding sites with improper balance of readpairs in double format. [2.0]
-samtools <arg> Path to samtools binary if not in PATH (i.e. /bin/samtools) [null].
-sr <arg> Filter sites with less than X SRs during breakpoint ascertainment. Default, 0, is to not filter any such sites. [0]
-t <arg> Path to the transposon ZIP file(s) to be used for this analysis.
-w <arg> Path to the working root directory.
-z <arg> Maximum reads to load into memory when iterating over sequence files. Setting higher increases run time, but may increase sensitivity in large (>60X coverage) bam files. Setting lower may decrease sensitivity in all bam files. [5000]
-help will print this message and exit
2017年６月現在、javaはversion1.８が最新になっている。ランするにはjava 1.7が必要になってくるので注意してください（java7リンク）。=> その後1.8でも動作はしました。
実行方法
入力はbwa memかalnで作成されたアライメントファイルである。 結果は標準的なvcf4.1のフォーマットで出力される。
マニュアルに記載されているヒトゲノムのデータでテストランしてみる。ラン前にhg19のbamとリファンレスゲノムをダウンロードする。
wget ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/phase3/data/NA12878/alignment/NA12878.mapped.ILLUMINA.bwa.CEU.low_coverage.20121211.bam
#bam.bai
wget ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/phase3/data/NA12878/alignment/NA12878.mapped.ILLUMINA.bwa.CEU.low_
#fasta
wget ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/technical/reference/phase2_reference_assembly_sequence/hs37d5.fa.gz
wget ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/technical/reference/phase2_reference_assembly_sequence/hs37d5.fa.gz.fai
#解凍
gzip -dv hs37d5.fa.gz
#リネーム
mv hs37d5.fa.gz.fai hs37d5.fa.fai
トランスポゾンの配列を準備する。テストデータにはALU、SVA、LINE1の配列が用意されている。その名前をlsコマンドで出力。
ls /full/path/to/MELTvX.X/me_refs/1KGP_Hg19/*.zip > mei_list.txt
exampleのトランスポゾン配列はzipで圧縮されている。
fastaで配列は記載され、bowtie2とsamtoolsでindexをつけられている。
準備ができたらランする。
java –jar MELT.jar Single -a –c 8 –h hs37d5.fa –l NA12878.mapped.ILLUMINA.bwa.CEU.low_coverage.20121211.bam –n add_bed_files/1KGP_Hg19/hg19.genes.bed –t mei_list.txt –w working_dir/
--single シングルゲノムの解析
perlで構築されているため、実行速度は遅めである。hg19の予測では、一晩走っていた。
指定したフォルダに結果は出力される。今回はALUS、SVA、LINE1の配列を指定したので、３つの配列の挿入部位がVCF4.1形式でコールされた。中身を確認する。
後半のフィールドが長くて見えづらい。前半のみ表示
SVA_final.vcf
LINE1_final.vcf
ALUは442コールあった。多すぎるので省略する。
トランスポソンタギングしたような株ではどうなるのかもテストしてみる。
検証中
注意点
トランスポゾン配列はfastaで指定する。名前はなんでもよいが、空白は許されていない。塩基配列は小文字、大文字とも認識可能で、60文字ごとに改行されていても１行に全部書かれていても動作する。ただしNがあってはいけない。MELTはATGCのみ認識する。また、トランスポゾン配列は１つずつファイルにする必要がある。マルチファスタ形式は許されていない。
repeatのBEDファイルが必要である。そのため、UCSCのcustom track機能で
frepeat masker(11) trackをつけてダウンロードする方法が推奨されている。ヒトゲノムならば、UCSCのtable browser（Table Browser）から、GROUPS => Repeats、track => RepeatMaskerにチェックをつけてダウンロードする。
引用
An integrated map of structural variation in 2,504 human genomes.
Nature. 2015 Oct 1;526(7571):75-81. doi: 10.1038/nature15394.
Sudmant PH1, Rausch T2, Gardner EJ3, Handsaker RE4,5, Abyzov A6, Huddleston J1,7, Zhang Y8,9, Ye K10,11, Jun G12,13, Hsi-Yang Fritz M2, Konkel MK14, Malhotra A15, Stütz AM2, Shi X16, Paolo Casale F17, Chen J8,18, Hormozdiari F1, Dayama G19, Chen K20, Malig M1, Chaisson MJ1, Walter K21, Meiers S2, Kashin S4,5, Garrison E22, Auton A23, Lam HY24, Jasmine Mu X8,25, Alkan C26, Antaki D27, Bae T6, Cerveira E15, Chines P28, Chong Z20, Clarke L17, Dal E26, Ding L10,11,29,30, Emery S31, Fan X20, Gujral M27, Kahveci F26, Kidd JM12,31, Kong Y23, Lameijer EW32, McCarthy S21, Flicek P17, Gibbs RA33, Marth G22, Mason CE34,35, Menelaou A36,37, Muzny DM38, Nelson BJ1, Noor A27, Parrish NF39, Pendleton M38, Quitadamo A16, Raeder B2, Schadt EE38, Romanovitch M15, Schlattl A2, Sebra R38, Shabalin AA40, Untergasser A2,41, Walker JA14, Wang M33, Yu F33, Zhang C15, Zhang J8,9, Zheng-Bradley X17, Zhou W20, Zichner T2, Sebat J27, Batzer MA14, McCarroll SA4,5; 1000 Genomes Project Consortium, Mills RE19,31, Gerstein MB8,9,42, Bashir A38, Stegle O17, Devine SE3, Lee C15,43, Eichler EE1,7, Korbel JO2,17.
https://www.ncbi.nlm.nih.gov/pubmed/26432246?dopt=Abstract
