遺伝性変異は複雑な遺伝性疾患における主要な寄与因子であると推測されている。多くの遺伝性疾患の遺伝率は比較的高いと推定されている。例えば、自閉症スペクトラム障害(ASD)の遺伝率は0.6を超え、統合失調症の遺伝率は0.5を超える。この遺伝性のごく一部のみが既知の遺伝性変異によって説明されており、"missing heritability"と名付けられた現象が存在する(ref.3 link)。一つの仮説はde novo突然変異である(特にindelと構造変異)。しかし、特にSV発見の複雑さから、これらの障害に対するそれらの寄与の説明は不完全になっている。一般的な遺伝性変異、特にde novoバリアントの発見は、依然として強い研究関心のあるトピックである。複雑な病因における遺伝性変異の役割を明らかにすることに加えて、多くのサンプルまたはコホートにわたるde novoバリアントの発見およびカタログ化の改善は、ヒトゲノムにおける重要な未解決問題にさらなる光を投げかける。
単一家系の全ゲノムシークエンシング(遺伝的疾患の孤立した症例を提示)は、生殖系列のde novo変異の発見から証明された成功例である。次に、ミスマッチ、ギャップ、カバレッジの急激なシフト、ペアリードのインサートサイズやペアリードの向きの不一致など、リードのアライメントに見られるアーティファクトに基づいて、各個人のバリアントが予測される。これは、継承されたバリエーション、真のデノボバリエーション、および偽のバリアントコールと識別するために検査する必要がある。
リファレンスに基づいた変異発見法は複雑な遺伝性疾患の研究において価値があることが証明されているが、本著者らはそれらの限界のいくつかに注目する。リードアラインメントアルゴリズムの一貫した改善にもかかわらず、各リードについて正しいマッピングを見つけることは、シーケンシングエラー、 repetitive DNA含有量、およびリファレンスにおけるミスアセンブリによって依然として複雑である。また、変異ブレイクポイントにまたがるか、または新規配列を含むためにリファレンスゲノムにマッピングされないリードは、マッピングベースの変異予測器では完全に無視される。最後に、リードアラインメントの分析によって決定されたほとんどのバリアントコールは、関心のある個人(子供または発端者)に固有のものではなく、代わりにファミリーとリファレンスゲノムドナーとの間の先祖の相違を反映している。ヒト生殖細胞系突然変異率の推定値は、一世代あたりおよそ80の新規突然変異という予想を与え、そして何百万もの遺伝的または偽のバリアントから真のde novo変異イベントを区別することは大きな課題である。
より一般的には、正確かつ包括的なde novoバリアントの発見は、いくつかの計算上および生物学的要因により複雑であり、依然としてとらえどころのない目標である。proband(発端者)にバリアントが存在することだけでなく、両親には存在しないことについてどのアルゴリズムも確信している必要がある。そして、一塩基変異(SNV)が最も一般的な変異のタイプであるが、よりまれなより大きな変異は全体的により多くのヌクレオチドに影響を及ぼし、そして遺伝性障害ではより大きな効果を有すると仮定される。 indelとSVの予測の本質的な複雑さのために、これらのより大きなde novoバリアントの正確な発見は特にチャレンジングになる。リファレンスマッピングのコンテキストでは、indelを確実にコールするには、すべてのギャップを一貫してアライメントした状態で、各リードをindelにまたがって正確にマッピングする必要がある。これは短いindelsに対してのみ可能であり、エラーやポジションギャップを起こしがちである。したがって、長さ10 bpを超えるindelsの予測は非常に困難であり、高い偽陽性率および高い偽陰性率を伴うことが証明されている。さらに、SVの予測は、リードデプス変化またはリードペアのシグネチャのような間接的なシグネチャによってのみ可能である。これらのシグネチャは非常にノイズが多くなる可能性があり、偽陽性および偽陽性の予測率が高くなる。結果として、それらの発生率を含むデノボSVのいくつかの基本的な特性は未知のままである。逆位を伴うduplicationなど、より複雑なタイプのde novo SVを予測する方法も存在しないことに注意することが重要である。
新規のバリアント予測に関する課題の多くは、リファレンスゲノムを介して間接的にではなく、直接関連する個体間のシーケンシング内容を比較するアプローチによって軽減することができる。このようなアプローチでは、リードの調整は不要である。また、offターゲットの共有/継承されたバリエーションには影響されない。マッピングフリーのアプローチが必要とするのは、アライメントアーチファクトに関して定義されていないバリエーションのシグネチャである。このマッピングフリーのコンテキストでは、k-mer(固定長kのシーケンス)の分析を介して、シーケンスの内容と量に関してバリアントを定義する。最近、ガン腫瘍における体細胞変異の正確な発見および全エキソームシーケンシングデータの完全な発見のために同様のアイデアが利用されている(ref.18, 19)。
直感的には、de novo生殖細胞系変異は、親ゲノムと比較してコーラーに新しい配列内容をもたらすはずである。最も単純な場合でさえ、十分に大きなkの値を考えると、突然変異をまたぐ k-merの全部ではないにしても、ほとんどの場合、単一のヌクレオチド置換はユニークであるべきである。ちなみに、これはindelやさまざまなタイプの構造的なバリエーションなど、他のクラスのバリエーションにも当てはまる。そして、probandゲノムは十分深くサンプリングされており、よってこれらの新規なk-merはシーケンシングデータ中に豊富に存在すると期待される。したがって、単一のモデルを使用して、SNVおよびより大きなバリアント(indel、構造変化)の両方を同時に検出することが可能であるはずである。
この直感に基づいて、本著者らは、新規のバリアント発見問題のマッピングフリーな定式化に基づく新しい方法であるKevlarを開発した。Kevlarは、両親からのリードには事実上欠けているが、proband/子供のデータのみに有意に高い存在量を持つ「興味深い」k-merを特定するために、k-merの存在量を調べる。これらのk-merは、デノボバリアントの潜在的存在の指標である。次に、ペアリード間で共有される新規のk-merに基づいてリードをグループ分けする。このそれぞれは推定バリアントを反映している。その後、Kevlarは標準的なアルゴリズムを使用して、リードの各セットをコンティグにアセンブルし、アセンブルされたコンティグをリファレンスゲノムにアライメントさせて予備的なバリアントコールを行う。最後に、Kevlarはバリアント予測をスコアリングし、誤ってコールされた遺伝性変異と真の de novo変異を区別するために新規な確率モデルを使用する。
論文では、シミュレートされたデータとリアルデータを使いこの新しい方法の有用性を実証する。 KevlarがSNVおよびshort indel発見のためのクラス最高のツールと同様の予測精度を達成すると同時に、より大きなイベントを高い精度で予測することを示す。ヌクレオチドレベルの精度でブレイクポイントを定義し、大規模なSVイベントを正確に予測するケブラーの能力を実証する。
Kevlarはオープンソースソフトウェアプロジェクトとして入手可能で、Python API、コマンドラインインターフェース、または標準のSnakemakeワークフローを介して呼び出すことができる。
The Kevlar pipeline. Preprintより転載。
Documentation
https://kevlar.readthedocs.io/en/latest/
Quick start guide
Quick start — kevlar 0.7rc1+1.ge55ba25 documentation
Tutorial
https://github.com/kevlar-dev/kevlar/tree/master/kevlar/workflows/mark-I
インストール
依存
- python=3
The kevlar package requires Python 3 and has several dependencies that are not in the standard Python libraries.
- the networkx package
- the pysam module
- the pandas library
- the scipy library
- the intervaltree library
- the khmer package
Also, kevlar requires the bwa command to be callable from your $PATH environmental variable.
本体 Github
mamba create -n Kevlar python=3.6 -y
conda activate Kevlar
pip3 install pysam networkx pandas scipy intervaltree git+https://github.com/dib-lab/khmer.git
pip3 install biokevlar
#snakemakeも必要
pip3 install snakemake
依存も含めて正しくインストールされているかテスト
pip3 install pytest
pytest --pyargs kevlar.tests
=========================================================================================================== test session starts ============================================================================================================
platform linux -- Python 3.6.7, pytest-4.4.0, py-1.8.0, pluggy-0.9.0
rootdir: /
collected 373 items
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_alac.py ........................... [ 7%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_assemble.py ............... [ 11%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_augment.py ..... [ 12%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_call.py ........................... [ 19%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_cigar.py ......... [ 22%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_cli.py .................... [ 27%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_count.py ........................... [ 34%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_dist.py .......... [ 37%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_evaluate.py . [ 37%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_filter.py ...... [ 39%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_gentrio.py ............................ [ 46%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_localize.py ............................ [ 54%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_mutablestring.py ..... [ 55%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_mutate.py ......... [ 58%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_novel.py ................ [ 62%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_partition.py ....... [ 64%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_progress.py .. [ 64%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_readgraph.py .. [ 65%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_readpair.py ......... [ 67%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_reference.py ..... [ 69%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_seqio.py .............. [ 72%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_simlike.py .................................... [ 82%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_sketch.py ................. [ 87%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_split.py .. [ 87%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_timer.py . [ 87%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_unband.py ... [ 88%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_varfilter.py .... [ 89%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_varmap.py .......................... [ 96%]
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_vcf.py ............ [100%]
============================================================================================================= warnings summary =============================================================================================================
usr/local/lib/python3.6/dist-packages/kevlar/tests/test_dist.py::test_tsv
/usr/local/lib/python3.6/dist-packages/kevlar/tests/test_dist.py:120: FutureWarning: read_table is deprecated, use read_csv instead.
data = pandas.read_table(distfile.name, sep='\t')
-- Docs: https://docs.pytest.org/en/latest/warnings.html
================================================================================================= 373 passed, 1 warnings in 70.72 seconds ==================================================================================================
ヘルプ
> kevlar count -h
# kevlar count -h
usage: kevlar count [-h] [-k K] [-c C] [-M MEM] [--max-fpr FPR] [--mask MSK]
[--count-masked] [--num-bands N] [--band I] [-t T]
counttable seqfile [seqfile ...]
Compute k-mer abundances for the provided sample. Supports k-mer banding:
see http://kevlar.readthedocs.io/en/latest/banding.html for more details.
positional arguments:
counttable name of the file to which the output (a k-mer count
table) will be written; the suffix ".counttable" will
be applied if the provided file name does not end in
".ct" or ".counttable"
seqfile input files in Fastq/Fasta format
optional arguments:
-h, --help show this help message and exit
-k K, --ksize K k-mer size; default is 31
-c C, --counter-size C
number of bits to allocate for counting each k-mer;
options are 1 (max count: 1), 4 (max count: 15), and 8
(max count: 255); default is 8
-M MEM, --memory MEM memory to allocate for the count table
--max-fpr FPR terminate if the estimated false positive rate for any
sample is higher than "FPR"; default is 0.2
--mask MSK counttable or nodetable of k-mers to ignore when
counting k-mers
--count-masked by default, when a mask is provided k-mers in the mask
are ignored; this setting inverts the behavior so that
only k-mers in the mask are counted
--num-bands N number of bands into which to divide the hashed k-mer
space
--band I a number between 1 and N (inclusive) indicating the
band to be processed
-t T, --threads T number of threads to use for file processing; default
is 1
Example::
kevlar count --memory 500M case1.ct case1-reads.fastq
Example::
kevlar count --ksize 25 --memory 12G --max-fpr 0.01 --threads 8 \
proband.counttable \
proband-R1.fq.gz proband-R2.fq.gz proband-unpaired.fq.gz
> kevlar novel -h
# kevlar novel -h
usage: kevlar novel --case F [F ...] [--case-counts F [F ...]]
[--control F [F ...]] [--control-counts F [F ...]] [-x X]
[-y Y] [-M MEM] [--max-fpr FPR] [--num-bands N] [--band I]
[-o FILE] [--save-case-counts CT [CT ...]]
[--save-ctrl-counts CT [CT ...]] [-h] [-k K]
[--abund-screen INT] [-t T] [--skip-until ID]
Identify "interesting" (potentially novel) k-mers and output the
corresponding reads. Here we define "interesting" k-mers as those which are
high abundance in each case sample and effectively absent (below some
specified abundance threshold) in each control sample.
Case/control config:
Specify input files, as well as thresholds for selecting "interesting"
k-mers. A single pass is made over input files for control samples (to
compute k-mer abundances), while two passes are made over input files for
case samples (to compute k-mer abundances, and then to identify
"interesting" k-mers). The k-mer abundance computing steps can be skipped
if pre-computed k-mer abunandances are provided using the "--case-counts"
and/or "--control-counts" settings. If "--control-counts" is declared, then
all "--control" flags are ignored. If "--case-counts" is declared,
FASTA/FASTQ files must still be provided with "--case" for selecting
"interesting" k-mers and reads.
--case F [F ...] one or more FASTA/FASTQ files containing reads from a
case sample; can be declared multiple times
corresponding to multiple case samples, see examples
below
--case-counts F [F ...]
counttable file(s) corresponding to each case sample;
if not provided, k-mer abundances will be computed
from FASTA/FASTQ input; only one counttable per
sample, see examples below
--control F [F ...] one or more FASTA/FASTQ files containing reads from a
control sample; can be declared multiple times
corresponding to multiple control samples, see
examples below
--control-counts F [F ...]
counttable file(s) corresponding to each control
sample; if not provided, k-mer abundances will be
computed from FASTA/FASTQ input; only one counttable
per sample, see examples below
-x X, --ctrl-max X k-mers with abund > X in any control sample are
uninteresting; default is X=1
-y Y, --case-min Y k-mers with abund < Y in any case sample are
uninteresting; default is Y=6
-M MEM, --memory MEM total memory allocated to k-mer abundance for each
sample; default is 1M; ignored when pre-computed k-mer
abundances are supplied via counttable
--max-fpr FPR terminate if the expected false positive rate for any
sample is higher than the specified FPR; default is
0.2
K-mer banding:
If memory is a limiting factor, it is possible to get a linear decrease in
memory consumption by running `kevlar novel` in "banded" mode. Splitting
the hashed k-mer space into N bands and only considering k-mers from one
band at a time reduces the memory consumption to approximately 1/N of the
total memory required. This implements a scatter/gather approach in which
`kevlar novel` is run N times, after the results are combined using
`kevlar filter`.
--num-bands N number of bands into which to divide the hashed k-mer
space
--band I a number between 1 and N (inclusive) indicating the
band to be processed
Output settings:
-o FILE, --out FILE file to which interesting reads will be written;
default is terminal (stdout)
--save-case-counts CT [CT ...]
save the computed k-mer counts for each case sample to
the specified count table file(s)
--save-ctrl-counts CT [CT ...]
save the computed k-mer counts for each control sample
to the specified count table file(s)
Miscellaneous settings:
-h, --help show this help message and exit
-k K, --ksize K k-mer size; default is 31
--abund-screen INT discard reads with any k-mers whose abundance is < INT
-t T, --threads T number of threads to use for file processing; default
is 1
--skip-until ID when re-running `kevlar novel`, skip all reads in the
case input until read with name `ID` is observed
Example::
kevlar novel --out novel-reads.augfastq --case proband-reads.fq.gz \
--control father-reads-r1.fq.gz father-reads-r2.fq.gz \
--control mother-reads.fq.gz
Example::
kevlar novel --out novel-reads.augfastq.gz \
--control-counts father.counttable mother.counttable \
--case-counts proband.counttable --case proband-reads.fastq \
--ctrl-max 0 --case-min 10 --ksize 27
Example::
kevlar novel --out output.augfastq \
--case proband1.fq --case proband2.fq \
--control control1a.fq control1b.fq \
--control control2a.fq control2b.fq \
--save-case-counts p1.ct p2.ct --save-ctrl-counts c1.ct c2.ct
> kevlar filter -h
# kevlar filter -h
usage: kevlar [-h] [-v] [-l F] [--tee] cmd ...
kevlar: error: argument cmd: invalid choice: 'filter\udcc2' (choose from 'count', 'novel', 'filter', 'augment', 'assemble', 'mutate', 'gentrio', 'partition', 'localize', 'call', 'alac', 'varfilter', 'simlike', 'split', 'dist', 'unband')
> kevlar partition -h
# kevlar partition -h
usage: kevlar partition [-h] [-s] [--min-abund X] [--max-abund Y] [--no-dedup]
[--gml FILE] [--split OUTPREFIX] [-o FILE]
infile
Construct a graph to group reads by shared interesting k-mers. Nodes in the
graph represent reads, and edges between a pair of nodes indicate that the two
corresponding reads have one or more interesting k-mers in common. Connected
components in the undirected graph correspond to distinct variants (or
variant-related breakpoints).
positional arguments:
infile input reads in augmented Fast[q|a] format
optional arguments:
-h, --help show this help message and exit
-s, --strict require perfect identity between overlapping reads for
inclusion in the same partition; by default, only a
shared interesting k-mer is required
--min-abund X ignore k-mers with abundance lower than X; default is 2
--max-abund Y ignore k-mers with abundance higher than Y; default is
200
--no-dedup skip step to remove duplicates
--gml FILE write read graph to .gml file
--split OUTPREFIX write each partition to a separate output file, each
with a filename like "OUTPREFIX.cc#.augfastq.gz"
-o FILE, --out FILE output file; default is terminal (stdout)
> kevlar call -h
# kevlar call -h
usage: kevlar call [-A A] [-B B] [-O O] [-E E] [--gen-mask FILE]
[--mask-mem MEM] [--mask-max-fpr FPR] [-h] [-d]
[--no-homopoly-filter] [--refr FILE] [-o FILE] [-k K]
queryseq targetseq
Align variant-related reads to the reference genome and call the variant from
the alignment.
Required inputs:
queryseq contigs assembled by "kevlar assemble"
targetseq region of reference genome identified by "kevlar
localize"
Alignment scoring:
-A A, --match A match score; default is 1
-B B, --mismatch B mismatch penalty; default is 2
-O O, --open O gap open penalty; default is 5
-E E, --extend E gap extension penalty; default is 0
Mask generation settings:
--gen-mask FILE generate a nodetable containing all k-mers that span
any variant call
--mask-mem MEM memory to allocate for the node table
--mask-max-fpr FPR terminate if the estimated false positive rate is
higher than "FPR"; default is 0.01
Miscellaneous settings:
-h, --help show this help message and exit
-d, --debug show debugging output
--no-homopoly-filter by default, short indels adjacent to homopolymers are
filtered out; use this flag to disable that filter
--refr FILE reference genome indexed for BWA search; if provided,
mates of interesting reads will be used to diambiguate
multi-mapping contigs
-o FILE, --out FILE output file; default is terminal (stdout)
-k K, --ksize K k-mer size; default is 31
> kevlar alac -h
# kevlar alac -h
usage: kevlar alac [-p ID] [--max-reads N] [-z Z] [-d D] [-x X]
[--include REGEX] [--exclude REGEX] [-A A] [-B B] [-O O]
[-E E] [--gen-mask FILE] [--mask-mem MEM]
[--mask-max-fpr FPR] [-h] [-o FILE] [-i I] [-k K] [-t T]
infile refr
Assemble partitioned reads, localize to the reference genome, align the
assembled contig to the reference, and call variants.
Required inputs:
infile partitioned reads in augmented Fastq format
refr reference genome in Fasta format (indexed for bwa
search)
Read assembly:
-p ID, --part-id ID only process partition "ID" in the input
--max-reads N do not attempt to assemble any partitions with more
than N reads (default: 10000)
Target extraction:
-z Z, --seed-size Z seed size; default is 51
-d D, --delta D retrieve the genomic interval from the reference by
extending beyond the span of all k-mer starting
positions by D bp
-x X, --max-diff X split and report multiple reference targets if the
distance between two seed matches is > X; by default,
X is set dynamically for each partition and is equal
to 3 times the length of the longest contig in the
partition; each resulting bin specifies a reference
target sequence to which assembled contigs will
subsequently be aligned
--include REGEX discard alignments to any chromosomes whose sequence
IDs do not match the given pattern
--exclude REGEX discard alignments to any chromosomes whose sequence
IDs match the given pattern
Alignment scoring:
-A A, --match A match score; default is 1
-B B, --mismatch B mismatch penalty; default is 2
-O O, --open O gap open penalty; default is 5
-E E, --extend E gap extension penalty; default is 0
Mask generation settings:
--gen-mask FILE generate a nodetable containing all k-mers that span
any variant call
--mask-mem MEM memory to allocate for the node table
--mask-max-fpr FPR terminate if the estimated false positive rate is
higher than "FPR"; default is 0.01
Miscellaneous settings:
-h, --help show this help message and exit
-o FILE, --out FILE output file; default is terminal (stdout)
-i I, --min-ikmers I do not report calls that a supported by fewer than `I`
interesting k-mers
-k K, --ksize K k-mer size; default is 31
-t T, --threads T process T partitions at a time using T threads
> kevlar unband -h
# kevlar unband -h
usage: kevlar unband [-h] [-n N] [-o FILE] infile [infile ...]
When kevlar is run in k-mer banding mode, the same read will typically
appear in multiple output files, annotated with a different set of
potentially novel k-mers in each case. This command will consolidate these
duplicated records across files into a single non-redundant set of reads
with the complete set of novel k-mer annotations.
positional arguments:
infile input files in augmented Fasta/Fastq format
optional arguments:
-h, --help show this help message and exit
-n N, --n-batches N number of batches into which records will be split;
default is 16; N temporary files are created and each
record from the input is written to a temporary file
based on its read name; then each batch is loaded into
memory and duplicate records are resolved
-o FILE, --out FILE output file; default is terminal (stdout)
> kevlar augment -h
# kevlar augment -h
usage: kevlar augment [-h] [-o FILE] augseqs seqs
Internally, kevlar annotates sequences with "interesting k-mers" and uses
"augmented" Fastq and Fasta formats. Processing sequences with third-part
tools usually requires discarding these annotations. This command is used to
augment/reaugment a set of sequences using annotations from an already
augmented sequence file.
positional arguments:
augseqs augmented sequence file
seqs sequences to annotate
optional arguments:
-h, --help show this help message and exit
-o FILE, --out FILE output file; default is terminal (stdout)
> kevlar mutate -h
# kevlar mutate -h
usage: kevlar mutate [-h] [-o FILE] mutations genome
Apply the specified mutations to the genome provided.
positional arguments:
mutations mutations file
genome genome to mutate
optional arguments:
-h, --help show this help message and exit
-o FILE, --out FILE output file; default is terminal (stdout)
テストラン
1、両親、子のfastq.gz、リファレンスゲノムをダウンロードする。
curl -L https://s3-us-west-1.amazonaws.com/noble-trios/neon-mother-reads.fq.gz -o mother.fq.gz
curl -L https://s3-us-west-1.amazonaws.com/noble-trios/neon-father-reads.fq.gz -o father.fq.gz
curl -L https://s3-us-west-1.amazonaws.com/noble-trios/neon-proband-reads.fq.gz -o proband.fq.gz
curl -L https://s3-us-west-1.amazonaws.com/noble-trios/neon-refr.fa.gz -o refr.fa.gz
2、リファレンスのindexing
bwa index refr.fa.gz
3、snakemakeのconfigの用意
curl -L https://s3-us-west-1.amazonaws.com/noble-trios/helium-config.json \
| sed "s:/home/user/Desktop:$(pwd):g" > helium-config.json
4、snakemakeのパイプラインの実行
#test run (dry run)
git clone https://github.com/kevlar-dev/kevlar.git
snakemake --configfile kevlar/kevlar/workflows/mark-I/config.json --snakefile kevlar/kevlar/workflows/mark-I/Snakefile --cores 64 --directory dir/ --printshellcmds --dryrun calls
#configファイルを書き換えてから実行する(ncbi/UniVec.faも必要)
snakemake \
--snakefile kevlar/workflows/mark-I/Snakefile \
--configfile helium-config.json --cores 4 --directory workdir -p calls
チュートリアルには使用する上での注意点もいくつか書かれています。確認してから使うようにして下さい。
kevlar/kevlar/workflows/mark-I at master · kevlar-dev/kevlar · GitHub
引用
Kevlar: a mapping-free framework for accurate discovery of de novo variants
Daniel S. Standage, C. Titus Brown, Fereydoun Hormozdi
bioRxiv preprint first posted online Feb. 13, 2019
Kevlar: A Mapping-Free Framework for Accurate Discovery of De Novo Variants
Daniel S. Standage, C. Titus Brown, Fereydoun Hormozdiar
iScience. 2019 Aug 30; 18: 28–36.