search
LoginSignup
1

More than 5 years have passed since last update.

posted at

updated at

PythonでPlink

PythonでPLINKのデータファイル(.bed, .fam, .map)などを利用したい。pyplinkというモジュールがある。

インストール

pipでインストールできる

shell
pip install pyplink

使い方

カレントディレクトリに foo.bed foo.fam foo.bim などの一連のファイルが存在するとすると

python3
from pyplink import PyPlink
pyp = PyPlink("foo")

これでpypというオブジェクトができる。これは.bed, .fam, .bimファイルを統合したオブジェクトである。各種メンバ関数でそれぞれの情報にアクセスできる。

python3
pyp.get_fam()
pyp.get_nb_samples()
pyp.get_bim()
pyp.get_nb_markers()
python3
markerNames = pyp.get_bim().iloc[:,5]

マーカー名を指定してジェノタイプを得る。acgtにすると塩基情報を得ることができる。

python3
pyp.get_geno_marker(markerNames[0])
pyp.get_acgt_geno_marker(markerNames[0])

イテレータとしてマーカーID, ジェノタイプを得ることもできる。

python3
markers = ["rs7092431", "rs9943770", "rs1578483"]
for marker_id, genoypes in pyp.iter_geno_marker(markers):
  print(marker_id)
  print(genotypes, end="\n\n")

サンプルスクリプト

23染色体上のマーカーについて男性サンプルのジェノタイプをすべて得る

python3
for marker_ID, genotypes in pyp.iter_geno_marker(y_markers):
    male_genotypes = genotypes[males]
    print("{:d} total genotypes".format(len(genotypes)))
    print("{:d} genotypes for {:,d} males ({} on chr{} and position {:,d})".format(
        len(male_genotypes),
        males.sum(),
        marker_ID,
        all_markers.loc[marker_ID, "chrom"],
        all_markers.loc[marker_ID, "pos"],
    ))
    break

指定したマーカーのMinor allele frequencyとジェノタイプを得る

python3
founders = (all_samples.father == "0") & (all_samples.mother == "0")
markers = ["rs7092431", "rs9943770", "rs1587483"]

for marker_ID, genotypes in pyp.iter_geno_marker(markers):
    valid_genotypes = genotypes[founders.values & (genotypes != -1)]
    maf = valid_genotypes.sum()/(len(valid_genotypes)*2)
    print(marker_ID, round(maf, 6), sep="\t")
    print(genotypes)

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
What you can do with signing up
1