ウェビナー参加者よりQ&Aに寄せられた質問

Q1.gitでダウンロードしてからjupyter notebookまでの作業についてもう一度教えてください。
conda activate pags
で事前に用意したpags環境に入って、
# コード取得
git clone https://github.com/genome-sci/python_bioinfo_2022
# ディレクトリに移動
cd python_bioinfo_2022

ターミナルで
jupyter notebook
実行

その後、
http://127.0.0.1:8888 または http://localhost:8888
をブラウザで開く。

場合によっては、画面に表示されたtokenつきのURLで開きます。
http://127.0.0.1:8888/?token=50387e…(以下省略)
or
http://localhost:8888/?token=50387e…(以下省略)

git コマンドが使えない方は下記リンクからzip形式でダウンロードできるので、展開し、その後展開してできたディレクトリに移動してください。
https://github.com/genome-sci/python_bioinfo_2022/archive/refs/heads/master.zip

Q2.事前準備でscikit-miscのインストールだけがうまくいきませんでした。エラー内容は以下の通りです。

PackagesNotFoundError: The following packages are not available from current channels:
(中略)
To search for alternate channels that may provide the conda package you’re looking for, navigate to
https://anaconda.org
and use the search bar at the top of the page.

上記のとおりに調べましたが、インストールコードはご案内いただいたものと同じconda install -c conda-forge scikit-misc しかありませんでした。
何が原因かわからなかったので、ご教授いただければ幸いです。

WindowsやM1 Macではインストールがうまくいかないケースがあるようです。
事前準備資料のトラブルシュートをご覧ください。
https://github.com/genome-sci/python_bioinfo_2022/blob/master/PREPARATION.md
Q3.谷澤先生の後半の練習問題の解答を見逃したので再掲していただけると助かります。
以下の通りです
i = 0
for record in SeqIO.parse(gbk_file_name, “genbank”):
for feature in record.features:
if feature.type != “CDS”:
continue # CDSでない場合には処理をスキップ
locus_tag = feature.qualifiers.get(“locus_tag”, [“”])[0]
### productを取得する処理を追加
product = feature.qualifiers.get(“product”, [“unknown_protein”])[0]
### 塩基配列を取得する処理を追加
seq = feature.location.extract(record.seq)
i += 1
### print で出力する処理を追加
# print(locus_tag, product, seq)
print(“>” + locus_tag + ” ” + product)
print(seq)
if i == 5: # 動作確認のため5件のみ出力したところで処理を停止
break
break # 動作確認のため最初の配列 (レコード) のみ処理するようにしている。この行を削除すれば全レコードを対象に処理される。
Q4.jupyternoteでIn[*]のようになり、実行ができなくなりました。
計算中を示します。n_dartsで大きすぎる数を指定すると結果が得られるまでに時間がかかります。
あまりにも時間がかかるようでしたら、jupyter notebookのメニューからKernel –> Restart を実行してください。
Q5.TPMが0.5以下の遺伝子をカットする場合もあるらしいのですが、カウント数のcut offはどのように決めるのでしょうか?
統計的な差と生物学的に意味がある差は別です。あまりにも少なすぎるリード数しかマップされない遺伝子は現実的に実験系として何も解釈ができません。そのため、シーケンスしたリード数やリファレンスにマップされたリード数等を元に最低カウント数を研究ごとに判断しているというのが現状です。
Q6.jupyter notebookを遺伝研等のサーバ上で動かすことはできますか?もしくはdesktopで起動させて、サーバ上のデータに直接アクセスできますか?
遺伝研スパコンを含むサーバー上で実行することは可能です。こちらをご参照ください。
https://sc.ddbj.nig.ac.jp/software/jupyter_notebook/
遺伝研スパコンの場合、jupyter notebookを動かしている計算ノードに直接接続することができないので、SSHポートフォワードという方法を使って接続していますが、通常のサーバーであればjupyter notebook起動時に表示されるURLを直接開けばOKです (firewall等でポートが開かれていることが前提です)
Q7.log2FCの計算で、微小な数を加えていましたが、微小な数でよいのでしょうか。割った時に、過大なFC値となりそうな気がしました。逆に1を加えるなどではだめなのでしょうか。
リードカウントに多大な影響を与えないような数を加えておいて、全体の計算でエラーを出さない、というくらいの意味なので、1でやってみてあまりおかしな値にならなければよいかと。
Q8.RNAseqのところで遺伝子長による正規化が重要になるのはどういった場合でしょうか?
転写産物の全長をシーケンスする場合、まったく同じコピー数の発現があっても、長さが長い遺伝子からは多くのシーケンスリードが観測されやすく、短い遺伝子からはリードが少なくなります。なので、遺伝子の長さでノーマライズしてmRNAのコピー数としての解釈ができるようにする必要があります。一方で、全長をシーケンスしない場合、たとえばmRNAの3’末端だけシーケンスする10X社のシングルセルシーケンスの場合などは、そのようなバイアスは生じないので、長さによる補正は必要ありません。
Q9.x軸,y軸以外の,黒い枠線を消すことはできますか?
一番簡単なのはseabornを使う方法です。
sns.despine()
をグラフ描画の最後につけると自動的に右と上の線が消えます。
Q10.import joypyのところで ModuleNotFoundError がでました。
何かインストールする必要がありますか?
conda install -c conda-forge joypy
でインストール可能です。(事前資料作成後に追加になりました)
Q11.No module named ‘matplotlib_venn’
とエラーが出ました。
conda install -c conda-forge matplotlib-venn
でインストール可能です。(事前資料作成後に追加になりました)
Q12.wt vs. mutantではなく、異なるがん種の遺伝子発現を比較したいときもTPMが使えるのでしょうか?
基本的には(数学的には)その通りです。
ただし、あまりに実験条件やRNA抽出プロトコルが異なる場合、そもそものRNA構成全体のパターンが異なることがあって、TPMの違いがDEGとは解釈できないケースがあります。実験条件がコントロールされた比較可能なサンプルかどうかを検討することが必要です。
そのあたりを議論した論文をご参照ください。
https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7373998/
Q13.今回の検定でカイ二乗検定を用いたのはなぜでしょうか?
例えば、あるトランスクリプトについてwt3個体とmutant3個体をT検定やマンホイットニーのU検定で比較するのはあまり良くない手法になりますか?
サンプルの数が3よりももっと増えればそれらの手法も良いと思いますが、N=3ではそれらの検定の検出力は小さく、ほとんどの場合で有意差は無しとなります。
Q14.scanpyで10Xのデータを扱う際にはTPMのような補正をする根拠がない とおっしゃっていたと思いますが、理由をもう一度教えていただけないでしょうか。
RNAの3’側の一部の領域のみ20-30 baseくらいシーケンスするため、bulk RNA-seqのようにRNA全長を呼んでいるわけでは無いので、遺伝子の長さの補正が不要なためです。
Q15.途中までseuratで解析しているデータを、一部だけscanpyで解析したいというときは、やはりanndataオブジェクトをはじめから作る必要があるでしょうか?
プロットのデザインの好みで違うプログラムを使いたいときや、seuratで整備されていない解析を一部取り入れたい場合に、簡単に互換するような方法はありますでしょうか。
Seuratでanndataに変換することもできますし、両方互換のフォーマットも一応あります。