Usage:
  prefetch [options] <path/SRA file | path/kart file> [<path/file> ...]
  prefetch [options] <SRA accession>
  prefetch [options] --list <kart_file>

prefetchの引数に、Run ID（SRAアクセッション番号）を指定して実行するのがメインの使い方となる。

以下に、Run ID、SRR17327096での実行結果を示した。

#SRR17327096のダウンロード
prefetch SRR17327096

#2022-01-04T14:52:02 prefetch.2.11.2: Current preference is set to retrieve SRA 
#Normalized Format files with full base quality scores.
#2022-01-04T14:52:03 prefetch.2.11.2: 1) Downloading 'SRR17327096'...
#2022-01-04T14:52:03 prefetch.2.11.2: SRA Normalized Format file is being retrieved, 
#if this is different from your preference, it may be due to current file availability.
#2022-01-04T14:52:03 prefetch.2.11.2:  Downloading via HTTPS...
#2022-01-04T14:52:10 prefetch.2.11.2:  HTTPS download succeed
#2022-01-04T14:52:10 prefetch.2.11.2:  'SRR17327096' is valid
#2022-01-04T14:52:10 prefetch.2.11.2: 1) 'SRR17327096' was downloaded successfully
#2022-01-04T14:52:10 prefetch.2.11.2: 'SRR17327096' has 0 unresolved dependencies

#ヘッド表示
head ./SRR17327096/SRR17327096.sra 
#NCBI.sra???Oqlockܥ?a$md֥?am"cur֥?a$4??md5֥?a$?)tblԥ?amSEQUENCE֥?am??col֥?am?!B`ALTREAD֥?am?#EX{??data֥?a$8$@#idx֥?a$((idx0֥?a$idx1֥?a$?
#md֥?am"cur֥?a$jmd5֥?a$T?QUALITY֥?am?#EX{??data֥?a$(?X?idx֥?a$x(idx0֥?a$idx1֥?a$0?idx2֥?a$p?md֥?am"cur֥?a$
#?md5֥?a$h?READ֥?am?#EX{??data֥?a$D?AZ?idx֥?a$(idx0֥?a$idx1֥?a$??idx2֥?a$T?md֥?am"cur֥?a$??md5֥?a$|?
#SPOT_GROUP֥?am?#EX{??data֥?a$xG?cidx֥?a$P(idx0֥?a$idx1֥?a$??idx2֥?a$??md֥?am"cur֥?a$,jmd5֥?a$@?md֥?#am"cur֥?a$??8fmd5֥?a$?)00??Y????"??c???W??a2fa47b25ad3fc0d8ffac2743996d730 *md/cur
#000dbe6a12f6b54989dfd016f340d3bb340 *md/cur
#000?tX߯tx??MD5CNTXT1234?#Eg?????ܺ?vT2?tX߯tx???tX߯tx??MD5CNTXT1234?#Eg?????ܺ?vT2?tX߯tx???tX߯tx??MD5CNTXT1234?#Eg?????ܺ?vT2?tX߯tx????
#row-len?schematypeINSDC:2na:packedversion 1;typedef B1 INSDC:2na:packed[2];alias INSDC:2na:packed INSDC:dna:2na;alias INSDC:2na:packed #NCBI:2na;06488d0d60f6f5ba40f1f5e0a9252b2ec *md/cur
#d1f0cf0771514cd3964222dab43dd0a1 *idx
#1c7d9a18aba301aa0419d34a0c853fd9 *idx1
#d41d8cd98f00b204e9800998ecf8427e *idx0
#32ce55e80d76d7db348d80134f20ff85 *idx2

実行の結果、ローカルディレクトリに「SRR17327096」フォルダが生成され、その中に「SRR17327096.sra」がダウンロードされる。「SRR17327096.sra」は、約 17.7 MBのバイナリファイルである。

prefetchの代表的なオプション

ここでは、prefetchでよく使われるオプションを表にまとめた。

オプション	概要
`-h`	すべてのオプション、一般的な使用法、およびバージョン情報を表示
`-V`	プログラムのバージョンを表示
`-f`	オブジェクトのダウンロードを強制
`-l`	kartファイルの内容をリストアップ
`-s`	ターゲットファイルサイズで kart ファイルの内容をリストアップ
`-N`	ダウンロードする最小ファイルサイズをKB単位で指定
`-X`	ダウンロードする最大ファイルサイズをKB単位で指定
`-a`	ascp プログラムと秘密鍵ファイル (asperaweb_id_dsa.openssh) へのパス
`-p`	ダウンロードの進行状況を表示
`--option-file`	ファイルから、オプションやパラメータを読み込む

複数のSRAファイルを一度にダウンロードする

NCBIのSequence Read Archiveでは、 Studyごとに一連のAccession Listが紐づけられている。例えば、Run ID、SRR17327096が含まれる StudyのすべてのRunデータを取得することにしてみる。

Accession Listを取得するために、まずはNCBIのSRR17327096サイトに行く。そこで、Studyの横にある、All runsをクリックする。

次のページで、DownloadのAccession Listをクリックすると、一連のRun IDが列挙された、SRR_Acc_List.txtがダウンロードされる。

このBioProject、PRJNA792067では、トータル 55 Run分のデータが含まれている。次の画像では、SRR_Acc_List.txtの一部を表示させた。

55 Run分はトライアルにしては、ちょっと多過ぎるので、上から3 lineだけを取り出して、別名のファイル(e.g. SRR_Acc_ListR.txt)を作成する。

この時の注意ポイントとして、 txtファイルの最終行に、必ず空行を入れておくこと。

prefetchコマンド実行の際には、 --option-fileにAccession Listを指定することで、複数のSRAファイルを一度にダウンロードできる。

以下に、SRR17327076、SRR17327077、SRR17327078を指定した時の実行結果を示す。

prefetch --option-file SRR_Acc_ListR.txt 

#2022-01-04T15:13:53 prefetch.2.11.2: Current preference is set 
#to retrieve SRA Normalized Format files with full base quality scores.
#2022-01-04T15:13:54 prefetch.2.11.2: 1) Downloading 'SRR17327076'...
#2022-01-04T15:13:54 prefetch.2.11.2: SRA Normalized Format file 
#is being retrieved, if this is different from your preference, 
#it may be due to current file availability.
#2022-01-04T15:13:54 prefetch.2.11.2:  Downloading via HTTPS...
#2022-01-04T15:14:01 prefetch.2.11.2:  HTTPS download succeed
#2022-01-04T15:14:01 prefetch.2.11.2:  'SRR17327076' is valid
#2022-01-04T15:14:01 prefetch.2.11.2: 1) 'SRR17327076' was downloaded successfully
#2022-01-04T15:14:01 prefetch.2.11.2: 'SRR17327076' has 0 unresolved dependencies

#2022-01-04T15:14:02 prefetch.2.11.2: Current preference is set to 
#retrieve SRA Normalized Format files with full base quality scores.
#2022-01-04T15:14:03 prefetch.2.11.2: 2) Downloading 'SRR17327077'...
#2022-01-04T15:14:03 prefetch.2.11.2: SRA Normalized Format file is 
#being retrieved, if this is different from your preference, it may 
#be due to current file availability.
#2022-01-04T15:14:03 prefetch.2.11.2:  Downloading via HTTPS...
#2022-01-04T15:14:14 prefetch.2.11.2:  HTTPS download succeed
#2022-01-04T15:14:14 prefetch.2.11.2:  'SRR17327077' is valid
#2022-01-04T15:14:14 prefetch.2.11.2: 2) 'SRR17327077' was downloaded successfully
#2022-01-04T15:14:14 prefetch.2.11.2: 'SRR17327077' has 0 unresolved dependencies

#2022-01-04T15:14:15 prefetch.2.11.2: Current preference is set to 
#retrieve SRA Normalized Format files with full base quality scores.
#2022-01-04T15:14:16 prefetch.2.11.2: 3) Downloading 'SRR17327078'...
#2022-01-04T15:14:16 prefetch.2.11.2: SRA Normalized Format file is 
#being retrieved, if this is different from your preference, it may 
#be due to current file availability.
#2022-01-04T15:14:16 prefetch.2.11.2:  Downloading via HTTPS...
#2022-01-04T15:14:26 prefetch.2.11.2:  HTTPS download succeed
#2022-01-04T15:14:26 prefetch.2.11.2:  'SRR17327078' is valid
#2022-01-04T15:14:26 prefetch.2.11.2: 3) 'SRR17327078' was downloaded successfully
#2022-01-04T15:14:26 prefetch.2.11.2: 'SRR17327078' has 0 unresolved dependencies

実行後に、以下のように、3つのフォルダが生成される。ただ、別々のフォルダに保存されるので、まとめ作業を後でしないといけないかも。

補足

シェルスクリプトで、複数のFASTA形式ファイル（+ gz圧縮）をダウンロードする。

簡単なシェルスクリプトを書いて、for-ループを実行することで、 SRR_Acc_ListR.txtを読み込んで、複数のFASTA形式ファイルを取得することも可能である。

複数のFASTA形式ファイル(1行あたり60塩基)を取得して、 gz圧縮する場合のシェルスクリプトの実行コードを以下に示す。

#複数のFASTA取得(1行あたり60塩基) + gz圧縮
while read line; do; 
  fastq-dump --gzip --fasta 60 $line; 
  done < ./SRR_Acc_ListR.txt
  
#Read 100959 spots for SRR17327076
#Written 100959 spots for SRR17327076
#Read 112855 spots for SRR17327077
#Written 112855 spots for SRR17327077
#Read 101553 spots for SRR17327078
#Written 101553 spots for SRR17327078

実行後に、以下のように、3つの.fasta.gzファイルが生成される。

ゲノム解析の関連記事

skume.net

参考資料

skume.net

www.server-memo.net

2022-01-03

R/ShortReadパッケージを使って、FASTQ形式ファイルを読み込む

Mac R言語 Bioconductor R - Bioconductor FASTQ形式 ShortRead ShortRead - readFastq SRA Toolkit fastq-dump シークエンスデータの取得シークエンスデータ解析 R - バイオインフォ FASTQファイル R - par() par - family日本語設定 R - ggplot2

はじめに
関連パッケージの読み込み
SRA Toolkit / fastq-dump コマンドを使って、FASTQファイルを取得する。
FASTQファイルの読み込み
リード長を取り出す方法
クオリティスコアを可視化する
サイクルごとにクオリティスコアを可視化する
まとめ
- ゲノム解析の関連記事
- 参考資料

はじめに

今回は、MacOSX環境で、SRA Toolkitのfastq-dumpを用いて取得した、 FASTQ形式ファイルをRに読み込む方法を扱う。

FASTQ形式ファイルの読み込みには、 Bioconductorパッケージの1つである、 ShortReadパッケージのreadFastq関数を用いる。

readFastq関数は、指定のディレクトリ内にあるFASTQ形式のファイルのうち、ファイル名がパターンマッチするものをすべて読み込む。あるいは、1つのファイル名を与えて、読み込むんこともできる。読み込み結果は、ファイル内の配列と品質スコアをコンパクトな内部表現として返す。

それでは、早速、FASTQファイルの読み込みを試してみる。

SRA Toolkit / fastq-dump コマンドを使って、FASTQファイルを取得する。

SRA Toolkitの設定は、以前の記事を参考のこと。

skume.net

ターミナルあるいはRStudioのTerminalタブにて、以下を実行する。

# fastq-dumpのパス確認
which fastq-dump

#SRR17327096データの取得(約16MB)
fastq-dump SRR17327096
#Read 110774 spots for SRR17327096
#Written 110774 spots for SRR17327096

#ファイルの確認
ls | grep SRR17327096      
#SRR17327096.fastq

FASTQファイルの読み込み

それでは、readFastq関数で、 FASTQファイルを読み込む。

#読み込み
fq <- ShortRead::readFastq("SRR17327096.fastq")

#表示
fq
#class: ShortReadQ
#length: 110774 reads; width: 488 cycles

#シークエンスをみる
reads <- ShortRead::sread(fq)
reads

リード長を取り出す方法

#リード長を取り出す
widths <- as.data.frame(reads@ranges@width)

#head表示
head(widths)
#  reads@ranges@width
#1                488
#2                488
#3                488
#4                488
#5                488
#6                488

#カウント
table(widths)
#widths
#   488 
#110774

クオリティスコアを可視化する

ここでは、FASTQファイルのクオリティスコアを抽出して、リードごとの平均値をグラフ化する。

# クオリティスコアを得る
quals <- Biostrings::quality(fq)

#表示
quals
#class: FastqQuality
#quality:
#BStringSet object of length 110774:
#         width seq
#     [1]   488 FFFFFFFFFFFFFFFFFFFFFFFFFFFFF...FFFFFFFFF:FFFFFF,FFFFFFFFFF:
#     [2]   488 FFFFFFFFFFFFFFFFFFFFFFFFFFFFF...FFFFFFFFFFFFFFFFF:FFFFFF:FFF
#     [3]   488 FFFFFFFFFFFFFFFFFFFFFFFFFFFFF...,FFFFFFFFFFFFFFFFFFF:FF,FFFF
#     [4]   488 FFFFFFFFFFFFFFFFFFFFFFFFFFFFF...FF:FFFF:FFFFF,FF::FFFFFFF,FF
#     [5]   488 FFFFFFFFFFFFFFFFFFFFFFF:FFFFF...FFFFFFFFF,FFFFFFF:F,FFFFFFFF
#     ...   ... ...
#[110770]   488 FFFFFFFFFFFFFFFFFFFFFFFFFF:FF...FFFF:FFF,FFFFFFFFFF:FFFFF:FF
#[110771]   488 FFFFFFFFFFFFFFFFFFFFFFFFFFFFF...FFFFFFFFFFFFFFFFFFFFFFFFFFFF
#[110772]   488 FFFFFFFFFFFFFFFFFFFFFFFFFFFFF...FFFFFFFFFFFFFFFFFFF:FFFF:FFF
#[110773]   488 FFFFFFFFFFFFFFFFFFFFFFFFF:FFF...F:FFFFFFFFFFFFFF:FFFFFFFFFF:
#[110774]   488 FFFFFFFFFFFFFFFFFFFFFFFFFFFFF...FFFFFFFFFFFFFFFFFFFFFFFFFFFF

#オブジェクトを強制的に行列クラスに変換する
scores <- as(quals, 'matrix')

#行ごとで平均値を計算
Avgscore <- data.frame(x=rowMeans(scores, na.rm = T))
str(Avgscore)
#'data.frame': 110774 obs. of  1 variable:
# $ x: num  36.4 36.8 36.3 36.5 35.8 ...

#可視化
ggplot(Avgscore) + 
  geom_histogram(aes(x=x))

サイクルごとにクオリティスコアを可視化する

#リストに変換する
Dat <- NULL
for(n in 1:ncol(scores)){
  #n <- 1
  print(n)
  a <- mean(scores[,n])
  b <- sd(scores[,n])
  cc <- data.frame(Ave=a, sd=b)
  Dat <- rbind(Dat, cc)
}

str(Dat)
#'data.frame': 488 obs. of  2 variables:
# $ Ave: num  36.5 36.6 36.2 36.6 36.3 ...
# $ sd : num  3.17 2.8 3.79 2.5 3.4 ...

#可視化
par(family="HiraKakuProN-W3", lwd=1.5, xpd=F, 
    mgp=c(2.25, 1, 0), mai=c(0.75,0.75, 0.5, 0.5))
plot(c(1, nrow(Dat)), c(0, 45), type = "n", xlab="Cycle", ylab="Q-Score")
points(Dat[,1], pch=21, col = "red", cex=0.25)
for(m in 1:ncol(scores)){segments(m, Dat[m,1]+Dat[m,2], m, Dat[m,1]-Dat[m,2], col = "red", lwd=0.1) }

Cycle数 v.s. Q-Score（赤丸はQ-Scoreの平均値で、平均値±SD（標準偏差）でドット表示）

まとめ

NCBIアーカイブからFASTQファイルのダウンロードから、ファイル読み込みまでの内容を扱った。

Bioconductorプロジェクトから関数群が提供されていることにも少し感動する。

ゲノム解析の関連記事

skume.net

参考資料

www.rdocumentation.org

bioconductor.org

kasperdanielhansen.github.io

github.com

2022-01-02

Wikipedia 英語記事「COVID-19」の日本語訳を公開してみた【その3: COVID感染症の病態生理】

COVID-19 SARS-CoV-2 Wikipedia 英語記事新型コロナウイルス生命科学

はじめに
Pathophysiology / 病態生理
過去の記事

はじめに

この記事は、2021年12月30日現在のWikipedia 英語記事「COVID-19」を日本語訳したものである。

en.wikipedia.org

この情報は英語ページにはありますが、まだ日本語ページは存在していません。結構しっかりした内容だったので、このwikipedia記事の日本語訳を作成することにしました。

コロナ関連の情報はすぐに古くなりますので、その都度、新しい情報を確認してください。

Pathophysiology / 病態生理

SARS-CoV-2ウイルス（COVID-19）は、生体内の様々な細胞やシステムに感染する。 COVID-19は、上気道（副鼻腔、鼻、喉）および下気道（気管および肺）に感染することがよく知られている*1。肺は、COVID-19の影響を最も受ける臓器である。なぜなら、ウイルスは、肺のII型肺胞細胞の表面に最も豊富に存在する、アンジオテンシン変換酵素2（ACE2）の受容体を介して、宿主細胞にアクセスする*2。ウイルスは、「スパイク」と呼ばれる、特殊な表面糖タンパク質を用いて、ACE2受容体とコネクトして、宿主細胞に侵入する*3。

COVID-19の病原性。（1A）COVID-19ウイルスは、エンドサイトーシス（endocytosis）あるいはACE2受容体に結合して膜融合を介して、上皮細胞に侵入し、ウイルスのRNAを細胞質内に放出する。（1B）ウイルスRNAは、細胞のシステムを利用して、ウイルスの非構造タンパク質や構造タンパク質を翻訳し、そのRNAを複製する。（1C）ウイルスの構造タンパク質S、E、Mは粗面小胞体（RER）で合わさる。（1D）その後、ウイルス構造とヌクレオキャプシドが小胞体ゴルジ体中間体（ERGIC）で集合体になる。（1E）ゴルジ小胞にパックされた、新しいウイルスは細胞膜と融合し、エキソサイトーシスを介して放出される。（2）COVID-19感染により炎症因子が誘導され、マクロファージや樹状細胞の活性化を引き起こす。（3）主要組織適合性複合体IおよびII（MHC IおよびII）を介した、COVID-19ウイルスの抗原提示は、体液性および細胞性免疫を刺激し、サイトカインおよび抗体の産生をもたらす。（4）COVID-19の重症例では、ウイルスは下気道に到達し、II型肺細胞に感染して、アポトーシスと表面活性物質（surfactant）の損失を引き起こす。マクロファージと好中球が流入して、サイトカインストーム*4を引き起こす。毛細血管の漏出により、肺胞水腫（alveolar edema）が生じる。ヒアルロン酸膜が形成される。これらの病理学的変化のすべては、肺胞のダメージや崩壊をもたらし、ガス交換を障害する。（出典元）

過去の記事

skume.net

*1:Harrison AG, Lin T, Wang P (December 2020). "Mechanisms of SARS-CoV-2 Transmission and Pathogenesis". Trends in Immunology. 41 (12): 1100–1115. doi:10.1016/j.it.2020.10.004. PMC 7556779. PMID 33132005.

*2:Verdecchia P, Cavallini C, Spanevello A, Angeli F (June 2020). "The pivotal link between ACE2 deficiency and SARS-CoV-2 infection". European Journal of Internal Medicine. 76: 14–20. doi:10.1016/j.ejim.2020.04.037. PMC 7167588. PMID 32336612.

*3:Letko M, Marzi A, Munster V (April 2020). "Functional assessment of cell entry and receptor usage for SARS-CoV-2 and other lineage B betacoronaviruses". Nature Microbiology. 5 (4): 562–569. doi:10.1038/s41564-020-0688-y. PMC 7095430. PMID 32094589.

*4:感染症などによって、サイトカイン（IL-1，IL-6，TNF-αなど）が異常に上昇して、その作用が全身に及ぶことで、好中球の活性化、血液凝固機構活性化、血管拡張などを全身に起こり、多臓器不全にまで進行する。この状態を、サイトカインストーム（cytokine storm）という。