自由研究3(厚労省の新型コロナの数値分析)

自由研究3(厚労省の新型コロナの数値分析)

目次

1 概要

  • 新型コロナの被害数値、厚労省から出されてるがPDFで処理しにくい
    • 厚労省からcsv形式の情報も出てるがPDFのまとめの数値とは異るデータ

2 関連文書

2.1 コロナ被害関係

2.1.1 厚労省

3 厚労省の過去のPDFデータへのリンク

4 厚労省のPDFファイルから4つの表を抜き出し、データベースに保存

4.1 説明動画


4.2 厚労省のPDFファイルから4つの表を抜き出し、データベースに保存の手順

4.2.1 必要なPythonパッケージのインストール

  1. 処理用の virtualenv を作成して、その環境に入る
    virtualenv corona_kourousyou
    cd corona_kourousyou
    source bin/activate
    
    • tikaをインストール(PDFファイルから文字を取り出すライブラリ)
    pip install tika
    pip install beautifulsoup4
    pip install nkf
    pip install pandas
    pip install matplotlib
    pip install ipython
    

4.2.2 過去の厚労省のPDFファイルをダウンロードして保存する

  • 今回表が4つになった後のPDFファイルのみを処理するので、令和2年12月2日以降のPDFファイルをダウンロード。私はwgetを用いて行なったが、好きなツールでOK
  • 令和3年9月1日までのPDFファイルをダウンロード
wget -m -l 1 https://www.mhlw.go.jp/content/10906000/000702069.pdf https://www.mhlw.go.jp/content/10906000/000704074.pdf https://www.mhlw.go.jp/content/10906000/000706508.pdf https://www.mhlw.go.jp/content/10906000/000710954.pdf https://www.mhlw.go.jp/content/10906000/000713230.pdf https://www.mhlw.go.jp/content/10906000/000716059.pdf https://www.mhlw.go.jp/content/10906000/000719997.pdf https://www.mhlw.go.jp/content/10906000/000724450.pdf https://www.mhlw.go.jp/content/10906000/000729529.pdf https://www.mhlw.go.jp/content/10906000/000734265.pdf https://www.mhlw.go.jp/content/10906000/000737675.pdf https://www.mhlw.go.jp/content/10906000/000741558.pdf https://www.mhlw.go.jp/content/10906000/000744930.pdf https://www.mhlw.go.jp/content/10906000/000748700.pdf https://www.mhlw.go.jp/content/10906000/000751942.pdf https://www.mhlw.go.jp/content/10906000/000755959.pdf https://www.mhlw.go.jp/content/10906000/000759133.pdf https://www.mhlw.go.jp/content/10906000/000764521.pdf https://www.mhlw.go.jp/content/10906000/000767330.pdf https://www.mhlw.go.jp/content/10906000/000769397.pdf https://www.mhlw.go.jp/content/10906000/000771872.pdf https://www.mhlw.go.jp/content/10906000/000775287.pdf https://www.mhlw.go.jp/content/10906000/000776153.pdf https://www.mhlw.go.jp/content/10906000/000779060.pdf https://www.mhlw.go.jp/content/10906000/000782185.pdf https://www.mhlw.go.jp/content/10906000/000785178.pdf https://www.mhlw.go.jp/content/10906000/000787852.pdf https://www.mhlw.go.jp/content/10906000/000790442.pdf https://www.mhlw.go.jp/content/10906000/000793909.pdf https://www.mhlw.go.jp/content/10906000/000796884.pdf https://www.mhlw.go.jp/content/10906000/000800096.pdf https://www.mhlw.go.jp/content/10906000/000803404.pdf https://www.mhlw.go.jp/content/10906000/000807083.pdf https://www.mhlw.go.jp/content/10906000/000809637.pdf https://www.mhlw.go.jp/content/10906000/000813216.pdf https://www.mhlw.go.jp/content/10906000/000816691.pdf https://www.mhlw.go.jp/content/10906000/000818427.pdf https://www.mhlw.go.jp/content/10906000/000820306.pdf https://www.mhlw.go.jp/content/10906000/000823800.pdf https://www.mhlw.go.jp/content/10906000/000826874.pdf

4.2.3 PDFファイルからデータを抽出して、データベースにデータを出力するPythonスクリプトの作成手順

  • pdfのあるディレクトリに移動 ( cd www.mhlw.go.jp/content/10906000/ )
  • 以下の内容で、pdffiltertika.pyを作成
from tika import parser
from bs4 import BeautifulSoup
import re
import sys
import os
import nkf
import pickle
import pandas as pd
import sqlite3
import datetime
#import warnings
#warnings.simplefilter("ignore")
fname=sys.argv[1]
bname,ext=os.path.splitext(fname)
pdf = parser.from_file(fname, xmlContent=True)
xhtml_data = BeautifulSoup(pdf['content'])

# データの整形
def myfilter(data):
  ans=[]
  tans=None
  redatetime=re.compile("令和([0-9]+)年([0-9]+)月([0-9]+)日([0-9]+)時.*")
  for i in data:
    for j in i.text.strip().split("\n"):
      k=re.sub("^[^0-9\.]+","",j.strip())
      if re.match("^[0-9\. ]+$",k):
	if len(k)!=0 and len(k.split())>2:
	  print(k)
	  if len(k.split())>8:
	    ans.append(k.split())
      if j.find("令和")!=-1:
	tmp=nkf.nkf('-Z1w' , j.replace(" ","")).decode('utf-8')
	m=redatetime.match(tmp)
	if m:
	  print(m)
	  print(m.group(0))
	  print(m.group(1))
	  print(m.group(2))
	  print(m.group(3))
	  print(m.group(4))
	  print(m.group())
	  print(tmp)
	  print(j)
	  tans=(m.group(1),m.group(2),m.group(3),m.group(4))
  return (tans,ans)

# 2ページ目の表の処理
t1,ans1=myfilter(xhtml_data.find_all('div', attrs={'class': 'page'})[1].find_all('p'))
# 3ページ目の表の処理
t2,ans2=myfilter(xhtml_data.find_all('div', attrs={'class': 'page'})[2].find_all('p'))

# 文字列を数値化
for i in range(6):
  ans1[i]=list(map(lambda x:int(x),ans1[i]))
for i in range(6,9):
  ans1[i]=list(map(lambda x:float(x),ans1[i]))
for i in range(1):
  ans2[i]=list(map(lambda x:float(x),ans2[i]))
for i in range(1,3):
  ans2[i]=list(map(lambda x:int(x),ans2[i]))

# 資料から抜き出した年月日を数値化
t1=list(map(lambda x:int(x),t1))
t2=list(map(lambda x:int(x),t2))

# 資料から抜き出した年月日をdatetime型に変更
tt=datetime.date(t1[0]+2021-3,t1[1],t1[2])
# 資料から抜き出した年月日を %Y-%m-%d 形式の文字列に変換(ex. 2021-09-01)
tts=datetime.date(t1[0]+2021-3,t1[1],t1[2]).strftime("%Y-%m-%d")

# 変数t1,t2,ans1,ans2をpickeファイルに書き出し
with open(bname+".pickle","wb") as f:
  pickle.dump((t1,t2,ans1,ans2),f)


# データベースファイルの初期化
with sqlite3.connect("covid19_data.db") as conn:
  # テーブルの生成
  conn.execute("create table if not exists pcr(t date, age00 real, age10 real , age20 real , age30 real , age40 real , age50 real , age60 real , age70 real , age80 real, age real)")
  conn.execute("create table if not exists pcr_m(t date, age00 real, age10 real , age20 real , age30 real , age40 real , age50 real , age60 real , age70 real , age80 real, age real)")
  conn.execute("create table if not exists pcr_f(t date, age00 real, age10 real , age20 real , age30 real , age40 real , age50 real , age60 real , age70 real , age80 real, age real)")
  conn.execute("create table if not exists death(t date, age00 real, age10 real , age20 real , age30 real , age40 real , age50 real , age60 real , age70 real , age80 real, age real)")
  conn.execute("create table if not exists death_m(t date, age00 real, age10 real , age20 real , age30 real , age40 real , age50 real , age60 real , age70 real , age80 real, age real)")
  conn.execute("create table if not exists death_f(t date, age00 real, age10 real , age20 real , age30 real , age40 real , age50 real , age60 real , age70 real , age80 real, age real)")
  conn.execute("create table if not exists mortality(t date, age00 real, age10 real , age20 real , age30 real , age40 real , age50 real , age60 real , age70 real , age80 real, age real)")
  conn.execute("create table if not exists mortality_m(t date, age00 real, age10 real , age20 real , age30 real , age40 real , age50 real , age60 real , age70 real , age80 real, age real)")
  conn.execute("create table if not exists mortality_f(t date, age00 real, age10 real , age20 real , age30 real , age40 real , age50 real , age60 real , age70 real , age80 real, age real)")
  conn.execute("create table if not exists severe_r(t date, age real, age00 real, age10 real , age20 real , age30 real , age40 real , age50 real , age60 real , age70 real , age80 real, nc real, investigate real, undisclosed real)")
  conn.execute("create table if not exists severe(t date, age integer, age00 integer, age10 integer , age20 integer , age30 integer , age40 integer , age50 integer , age60 integer , age70 integer , age80 integer, nc integer, investigate integer, undisclosed integer)")
  conn.execute("create table if not exists hospitalization(t date, age integer, age00 integer, age10 integer , age20 integer , age30 integer , age40 integer , age50 integer , age60 integer , age70 integer , age80 integer, nc integer, investigate integer, undisclosed integer)")

# 取り込んだデータを書き出し
with sqlite3.connect("covid19_data.db") as conn:
  conn.execute("insert into pcr values( ?,?,?,?,? ,?,?,?,?,?, ? )", [tts]+ans1[0])
  conn.execute("insert into pcr_m values( ?,?,?,?,? ,?,?,?,?,?, ? )", [tts]+ans1[1])
  conn.execute("insert into pcr_f values( ?,?,?,?,? ,?,?,?,?,?, ? )", [tts]+ans1[2])
  conn.execute("insert into death values( ?,?,?,?,? ,?,?,?,?,?, ? )", [tts]+ans1[3])
  conn.execute("insert into death_m values( ?,?,?,?,? ,?,?,?,?,?, ? )", [tts]+ans1[4])
  conn.execute("insert into death_f values( ?,?,?,?,? ,?,?,?,?,?, ? )", [tts]+ans1[5])
  conn.execute("insert into mortality values( ?,?,?,?,? ,?,?,?,?,?, ? )", [tts]+ans1[6])
  conn.execute("insert into mortality_m values( ?,?,?,?,? ,?,?,?,?,?, ? )", [tts]+ans1[7])
  conn.execute("insert into mortality_f values( ?,?,?,?,? ,?,?,?,?,?, ? )", [tts]+ans1[8])
  conn.execute("insert into severe_r values( ?,?,?,?,? ,?,?,?,?,?, ?,?,?,? )", [tts]+ans2[0])
  conn.execute("insert into severe values( ?,?,?,?,? ,?,?,?,?,?, ?,?,?,? )", [tts]+ans2[1])
  conn.execute("insert into hospitalization values( ?,?,?,?,? ,?,?,?,?,?, ?,?,?,? )", [tts]+ans2[2])

# SQLでテスト用
# select * from death order by t asc;
  • 8個のスペースがtabに置き換わってしまってるので以下で変換(お好きなエディタで置換すればOK)
sed -e "s/\t/        /g" pdffiltertika.py > 1.py
mv 1.py pdffiltertika.py

4.3 作業自動化のためのGNU MakeのMakefileを作成およびデータ抽出

4.3.1 作業自動化のためのGNU MakeのMakefileを作成 手順

  • pdfのあるディレクトリに移動 ( cd www.mhlw.go.jp/content/10906000/ )
  • 以下の内容でMakefile作成
PDF=${wildcard *.pdf}
PICKLE=${PDF:.pdf=.pickle}

all: ${PICKLE}

%.pickle: %.pdf
	python3 pdffiltertika.py $<

4.3.2 データ抽出

  • データベースファイルの削除 と pickleファイルを削除
rm covid19_data.db *.pickle
  • PDFファイルからデータ抽出
make 

make

4.4 この章のまとめ

  • 厚労省のPDFファイルから、表の数値をプログラムで取り出せるようにした。

5 厚労省のPDFファイルから4つの表を抜き出し、データベースに保存したデータのpandasのデータフレームへ取り込み

5.1 説明動画


5.2 厚労省のPDFファイルから4つの表を抜き出し、データベースに保存したデータのpandasのデータフレームへ取り込み手順

5.2.1 Pythonの対話環境起動

  • virtualenv 環境にはいる
source bin/activate
cd  www.mhlw.go.jp/content/10906000/
  • Pythonの対話環境起動
ipython

5.2.2 SQLite3のデータベースファイルからデータ取り出し

  • 以下の名前の変数に時系列データを入れる
    • pcrに陽性数
      • pcr_mに陽性数(男性
      • pcr_fに陽性数(女性
    • deathに死者数の累計
      • death_mに死者数の累計(男性
      • death_fに死者数の累計(女性
    • mortalityに死亡
      • mortality_mに死亡率(男性
      • mortality_fに死亡率(女性
    • severe に重傷者数(累計ではないその時点の)
    • severe_r に重傷化率
    • hospitalization に入院を必要とする人数
import pandas as pd
import sqlite3
import datetime
sql = "select * from death order by t asc"
with sqlite3.connect("covid19_data.db") as conn:
  death=pd.read_sql(sql, conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_m=pd.read_sql(sql.replace("death","death_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_f=pd.read_sql(sql.replace("death","death_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality=pd.read_sql(sql.replace("death","mortality"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_m=pd.read_sql(sql.replace("death","mortality_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_f=pd.read_sql(sql.replace("death","mortality_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr=pd.read_sql(sql.replace("death","pcr"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_m=pd.read_sql(sql.replace("death","pcr_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_f=pd.read_sql(sql.replace("death","pcr_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe_r=pd.read_sql(sql.replace("death","severe_r"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe=pd.read_sql(sql.replace("death","severe"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  hospitalization=pd.read_sql(sql.replace("death","hospitalization"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})

5.2.3 取り込んだデータの確認

#死亡者数の累計表示
death
#1週間毎の死亡者数
death.diff()
#1週間毎の陽性者数あたりの死亡率
death.diff()/pcr.diff()*100.0

5.3 この章のまとめ

  • データベースに格納してあった、時系列データをPythonのpandasのデータフレームに取り出した

6 厚労省のPDFファイルから4つの表を抜き出したデータを色々みてみる(2021/9/1まで)

6.1 説明動画


6.2 厚労省のPDFファイルから4つの表を抜き出したデータを色々いじってみる手順

6.2.1 Pythonの対話環境起動

  • virtualenv 環境にはいる
source bin/activate
cd  www.mhlw.go.jp/content/10906000/
  • Pythonの対話環境起動
ipython

6.2.2 データ取り込み

import pandas as pd
import sqlite3
import datetime
sql = "select * from death order by t asc"
with sqlite3.connect("covid19_data.db") as conn:
  death=pd.read_sql(sql, conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_m=pd.read_sql(sql.replace("death","death_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_f=pd.read_sql(sql.replace("death","death_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality=pd.read_sql(sql.replace("death","mortality"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_m=pd.read_sql(sql.replace("death","mortality_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_f=pd.read_sql(sql.replace("death","mortality_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr=pd.read_sql(sql.replace("death","pcr"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_m=pd.read_sql(sql.replace("death","pcr_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_f=pd.read_sql(sql.replace("death","pcr_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe_r=pd.read_sql(sql.replace("death","severe_r"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe=pd.read_sql(sql.replace("death","severe"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  hospitalization=pd.read_sql(sql.replace("death","hospitalization"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})

6.2.3 最近多い重傷者の変化を調べてみる

  1. 小数点以下の表示桁数を変更
    pd.options.display.precision=3
    
  2. 重傷者数の時系列表示
    severe
    
  3. 重傷者数の変化時系列表示
    1. 重傷者数の1週間の変化を表示
      severe.diff()
      
  4. 重症化率を調べてみる
    severe_r
    
  5. 重傷化率の変化時系列
    1. 重傷者化率の1週間の変化を
      severe_r.diff()
      
  6. 死者数
    1. 死者数の時系列
      death
      
      • 最後の行を表示
      death.tail(1)
      
      • 累積和を計算
      death.tail(1).cumsum(axis=1)
      
      • 累積和の割合を計算
      death.tail(1).cumsum(axis=1)/float(death.tail(1).cumsum(axis=1)["age80"])*100
      
    2. 1週間毎の死者数の時系列
      death.diff()
      
    3. 1週間毎の死者数の4週移動平均
      death.diff().rolling(4).mean()
      
  7. 死者数の変化
    1. 1週間毎の死者数の変化の時系列表示
      death.diff().diff()
      
  8. 死亡率を調べてみる
    1. 1週間毎の死亡率を調べてみる
      death.diff()/pcr.diff()*100
      
    2. 1週間毎の死亡率の1ヶ月の移動平均
      (death.diff()/pcr.diff()*100).rolling(4).mean()
      
  9. 死亡率の変化を調べてみる
    1. 1週間毎の死亡率の変化を調べてみる
      (death.diff()/pcr.diff()*100).diff()
      
    2. 1週間毎の死亡率の1ヶ月の移動平均の差分
      (death.diff()/pcr.diff()*100).rolling(4).mean().diff()
      

7 Matplotlibでグラフ表示1(Matplotlibの基礎確認)

7.1 説明動画


7.2 matplotlibでグラフ表示の手順

7.2.1 Pythonの対話環境起動

  • virtualenv 環境にはいる
source bin/activate
cd  www.mhlw.go.jp/content/10906000/
  • matplotlibがはいってないなら以下でインストール
pip install matplotlib
  • Pythonの対話環境起動
ipython

7.2.2 データ取り込み等

import pandas as pd
import sqlite3
import datetime
import matplotlib.pyplot as plt
import matplotlib
sql = "select * from death order by t asc"
with sqlite3.connect("covid19_data.db") as conn:
  death=pd.read_sql(sql, conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_m=pd.read_sql(sql.replace("death","death_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_f=pd.read_sql(sql.replace("death","death_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality=pd.read_sql(sql.replace("death","mortality"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_m=pd.read_sql(sql.replace("death","mortality_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_f=pd.read_sql(sql.replace("death","mortality_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr=pd.read_sql(sql.replace("death","pcr"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_m=pd.read_sql(sql.replace("death","pcr_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_f=pd.read_sql(sql.replace("death","pcr_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe_r=pd.read_sql(sql.replace("death","severe_r"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe=pd.read_sql(sql.replace("death","severe"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  hospitalization=pd.read_sql(sql.replace("death","hospitalization"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})

7.2.3 重傷者の数をグラフ表示1

severe.plot()
plt.show()

7.2.4 重傷者の数をグラフ表示2 グリッドを追加して表示

severe.plot()
plt.grid()
plt.show()

7.2.5 日本語フォントを指定

  • これしないと文字化けする
  • フォントのフルパスは実際使いたいフォントパスに変更してください。(“/usr/share/fonts/truetype/fonts-japanese-gothic.ttf”は私のUbuntu環境のfontのパスから適当に選択しました)
matplotlib.rcParams["font.family"] = matplotlib.font_manager.FontProperties(fname="/usr/share/fonts/truetype/fonts-japanese-gothic.ttf").get_name()

7.2.6 重傷者の数をグラフ表示3 軸ラベルも追加して表示

severe.plot()
plt.xlabel("年月日")
plt.ylabel("重傷者数(人)")
plt.title("厚労省のPDFから抜き出した重傷者数(人)")
plt.grid()
plt.show()

7.3 この章のまとめ

  • Matplotlibの使い方の基礎確認

8 Matplotlibでグラフ表示2(上下に2つのグラフ表示)

8.1 説明動画


8.2 matplotlibでグラフ表示の手順

8.2.1 Pythonの対話環境起動

  • virtualenv 環境にはいる
source bin/activate
cd  www.mhlw.go.jp/content/10906000/
  • matplotlibがはいってないなら以下でインストール
pip install matplotlib
  • Pythonの対話環境起動
ipython

8.2.2 データ取り込み等

import pandas as pd
import sqlite3
import datetime
import matplotlib.pyplot as plt
import matplotlib
sql = "select * from death order by t asc"
with sqlite3.connect("covid19_data.db") as conn:
  death=pd.read_sql(sql, conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_m=pd.read_sql(sql.replace("death","death_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_f=pd.read_sql(sql.replace("death","death_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality=pd.read_sql(sql.replace("death","mortality"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_m=pd.read_sql(sql.replace("death","mortality_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_f=pd.read_sql(sql.replace("death","mortality_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr=pd.read_sql(sql.replace("death","pcr"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_m=pd.read_sql(sql.replace("death","pcr_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_f=pd.read_sql(sql.replace("death","pcr_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe_r=pd.read_sql(sql.replace("death","severe_r"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe=pd.read_sql(sql.replace("death","severe"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  hospitalization=pd.read_sql(sql.replace("death","hospitalization"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})

8.2.3 重傷者の数をグラフ表示3 軸ラベルも追加して表示

fig = plt.figure()
axes= fig.subplots(2)
severe.plot(ax=axes[0])
death.diff().plot(ax=axes[1])
axes[0].grid()
axes[0].set_xlabel("年月日")
axes[0].set_ylabel("重傷者数(人)")
axes[0].set_title("厚労省PDFから抜き出した数値から作成")
axes[1].grid()
axes[1].set_xlabel("年月日")
axes[1].set_ylabel("1週間毎の死者数(人)")
plt.show()

8.3 この章のまとめ

  • Matplotlibの使い方の基礎、複数のグラフの同時表示のやりかたの確認をした

9 Matplotlibでグラフ表示3(パイチャート)

9.1 説明動画


9.2 matplotlibでグラフ表示の手順

9.2.1 Pythonの対話環境起動

  • virtualenv 環境にはいる
source bin/activate
cd  www.mhlw.go.jp/content/10906000/
  • matplotlibがはいってないなら以下でインストール
pip install matplotlib
  • Pythonの対話環境起動
ipython

9.2.2 データ取り込み等

  • フォントのフルパスは実際使いたいフォントパスに変更してください。(“/usr/share/fonts/truetype/fonts-japanese-gothic.ttf”は私のUbuntu環境のfontのパスから適当に選択しました)
import pandas as pd
import sqlite3
import datetime
import matplotlib.pyplot as plt
import matplotlib
sql = "select * from death order by t asc"
with sqlite3.connect("covid19_data.db") as conn:
  death=pd.read_sql(sql, conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_m=pd.read_sql(sql.replace("death","death_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_f=pd.read_sql(sql.replace("death","death_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality=pd.read_sql(sql.replace("death","mortality"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_m=pd.read_sql(sql.replace("death","mortality_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_f=pd.read_sql(sql.replace("death","mortality_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr=pd.read_sql(sql.replace("death","pcr"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_m=pd.read_sql(sql.replace("death","pcr_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_f=pd.read_sql(sql.replace("death","pcr_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe_r=pd.read_sql(sql.replace("death","severe_r"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe=pd.read_sql(sql.replace("death","severe"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  hospitalization=pd.read_sql(sql.replace("death","hospitalization"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
matplotlib.rcParams["font.family"] = matplotlib.font_manager.FontProperties(fname="/usr/share/fonts/truetype/fonts-japanese-gothic.ttf").get_name()

9.2.3 パイチャート

  1. 死者数の累計
    tmp=death.tail(1).T[:-1]
    tmp.index=["10歳未満","10代","20代","30代","40代","50代","60代","70代","80代以上"]
    tmp.plot(kind='pie', autopct='%1.1f%%',subplots=True)
    plt.title("新型コロナでの死者累計の割合2021/09/01")
    plt.show()
    
  2. 直近重傷者
    tmp2=severe.tail(1).T[1:-3]
    tmp2.index=["10歳未満","10代","20代","30代","40代","50代","60代","70代","80代以上"]
    tmp2.plot(kind='pie', autopct='%1.1f%%',subplots=True)
    plt.title("新型コロナでの重傷者の割合2021/09/01")
    plt.show()
    
  3. 直近1週間の死者
    tmp3=death.diff().tail(1).T[:-1]
    tmp3.index=["10歳未満","10代","20代","30代","40代","50代","60代","70代","80代以上"]
    tmp3.plot(kind='pie', autopct='%1.1f%%',subplots=True)
    plt.title("新型コロナでの直近1週間の死者の割合2021/09/01")
    plt.show()
    
  4. 3つのグラフ同時表示
    fig = plt.figure()
    axes= fig.subplots(1,3)
    tmp.plot(kind='pie',ax=axes[0], autopct='%1.1f%%',subplots=True,legend=True)
    tmp2.plot(kind='pie',ax=axes[1], autopct='%1.1f%%',subplots=True,legend=False)
    tmp3.plot(kind='pie',ax=axes[2], autopct='%1.1f%%',subplots=True,legend=False)
    axes[0].set_title("新型コロナでの死者累計の割合2021/09/01")
    axes[1].set_title("新型コロナでの重傷者の割合2021/09/01")
    axes[2].set_title("新型コロナでの直近1週間の死者の割合2021/09/01")
    plt.show()
    
    1. 東京都福祉保健局(基礎疾患の有無、被害者の平均年齢情報あり)

9.3 この章のまとめ

  • Matplotlibの使い方の基礎、パイチャートの描画

10 新規PDFデータ追加

10.1 説明動画


10.2 データ追加の手順

10.2.1 virtualenv 環境にはいる

source bin/activate

10.2.2 pdfファイル追加

  • 以下は令和3年9月8日のPDFファイルのURLなので、更新された分(というか処理したいPDFを入れて)処理する必要がある
  • 複数ファイル追加もOK
  • ただし、令和2年12月2日より前のPDFファイルはフォーマットが異るため、入れるとエラーになる
  • 今後フォーマットが変更されると同じ処理で上手くデータ追加出来なくなる
wget -m -l 1 https://www.mhlw.go.jp/content/10906000/000830345.pdf

10.2.3 データ追加

  • pdfファイルのあるディレクトリに移動
cd  www.mhlw.go.jp/content/10906000/
  • GNU makeで処理
    • これで新しく追加した分のPDFファイルを処理して、抜き出した情報をデータベースに格納
make

10.3 この章のまとめ

  • 新しくPDFファイルを追加したときの処理方法についての説明

11 主要なデータの確認及びグラフ表示

11.1 説明動画


11.2 主要なデータの確認及びグラフ表示

11.2.1 データ取り込み等

  • フォントのフルパスは実際使いたいフォントパスに変更してください。(“/usr/share/fonts/truetype/fonts-japanese-gothic.ttf”は私のUbuntu環境のfontのパスから適当に選択しました)
  • 以下の名前の変数に時系列データを入れる
    • pcrに陽性数
      • pcr_mに陽性数(男性
      • pcr_fに陽性数(女性
    • deathに死者数の累計
      • death_mに死者数の累計(男性
      • death_fに死者数の累計(女性
    • mortalityに死亡
      • mortality_mに死亡率(男性
      • mortality_fに死亡率(女性
    • severe に重傷者数(累計ではないその時点の)
    • severe_r に重傷化率
    • hospitalization に入院を必要とする人数
import pandas as pd
import sqlite3
import datetime
import matplotlib.pyplot as plt
import matplotlib
sql = "select * from death order by t asc"
with sqlite3.connect("covid19_data.db") as conn:
  death=pd.read_sql(sql, conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_m=pd.read_sql(sql.replace("death","death_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  death_f=pd.read_sql(sql.replace("death","death_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality=pd.read_sql(sql.replace("death","mortality"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_m=pd.read_sql(sql.replace("death","mortality_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  mortality_f=pd.read_sql(sql.replace("death","mortality_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr=pd.read_sql(sql.replace("death","pcr"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_m=pd.read_sql(sql.replace("death","pcr_m"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  pcr_f=pd.read_sql(sql.replace("death","pcr_f"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe_r=pd.read_sql(sql.replace("death","severe_r"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  severe=pd.read_sql(sql.replace("death","severe"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
  hospitalization=pd.read_sql(sql.replace("death","hospitalization"), conn, index_col='t', parse_dates={'Date': '%Y-%m-%d'})
# 日本語フォントの設定
matplotlib.rcParams["font.family"] = matplotlib.font_manager.FontProperties(fname="/usr/share/fonts/truetype/fonts-japanese-gothic.ttf").get_name()
# 小数点以下の表示桁数を変更
pd.options.display.precision=3

11.2.2 死者数、最新から1週間での死者数、重傷者のパイチャート表示

fig = plt.figure()
axes= fig.subplots(1,5)
tmp=death.tail(1).T[:-1]
tmp2=severe.tail(1).T[1:-3]
tmp3=death.diff().tail(1).T[:-1]
tmp4=hospitalization.tail(1).T[1:-3]
tmp.index=["10歳未満","10代","20代","30代","40代","50代","60代","70代","80代以上"]
tmp2.index=["10歳未満","10代","20代","30代","40代","50代","60代","70代","80代以上"]
tmp3.index=["10歳未満","10代","20代","30代","40代","50代","60代","70代","80代以上"]
tmp4.index=["10歳未満","10代","20代","30代","40代","50代","60代","70代","80代以上"]
tmp.plot(kind='pie',ax=axes[0], autopct='%1.1f%%',subplots=True,legend=False)
tmp2.plot(kind='pie',ax=axes[1], autopct='%1.1f%%',subplots=True,legend=False)
tmp3.plot(kind='pie',ax=axes[2], autopct='%1.1f%%',subplots=True,legend=False)
tmp4.plot(kind='pie',ax=axes[3], autopct='%1.1f%%',subplots=True,legend=False)
tmp.plot(kind='pie',ax=axes[4], autopct='%1.1f%%',subplots=True,legend=True)
axes[0].set_title("死者累計の割合")
axes[1].set_title("重傷者の割合")
axes[2].set_title("直近1週間の死者の割合")
axes[3].set_title("入院治療等を要する者")
plt.show()

11.2.3 1週間の死者数の計算およびグラフ表示

  1. 重傷者数
    severe
    
  2. 1週間の死者数
    death.diff()
    
  3. グラフ表示
    fig = plt.figure()
    axes= fig.subplots(2)
    severe.plot(ax=axes[0])
    death.diff().plot(ax=axes[1])
    axes[0].grid()
    axes[0].set_xlabel("年月日")
    axes[0].set_ylabel("重傷者数(人)")
    axes[0].set_title("厚労省PDFから抜き出した数値から作成")
    axes[1].grid()
    axes[1].set_xlabel("年月日")
    axes[1].set_ylabel("1週間毎の死者数(人)")
    plt.show()
    

11.2.4 1週間毎の死亡率

  1. 計算
    death.diff()/pcr.diff()*100
    
  2. グラフ表示
    (death.diff()/pcr.diff()*100).plot()
    plt.title("1週間毎の死者数(人)")
    plt.show()
    

11.2.5 死亡率(その時までの累計での)

  1. データ確認
    mortality
    
  2. グラフ表示
    mortality.plot()
    plt.title("死者率(%)")
    plt.show()
    

11.2.6 重症化率

  1. データ確認
    severe_r
    
  2. グラフ表示
    severe_r.plot()
    plt.title("重症化率(%)")
    plt.show()
    

11.2.7 入院治療等を要する者(人数)

  1. データ確認
    hospitalization
    
  2. グラフ表示
    hospitalization.plot()
    plt.title("入院治療等を要する者(人)")
    plt.show()
    

11.3 この章のまとめ

  • 主要な情報の確認と、グラフ化をおこなった

12 今後

  • 今後文書追加するかも、しないかも

13 この文書のチェンジログ

  • 2021/09/02 初版
  • 2021/09/06 修正
  • 2021/09/08 以下の章を追加
    • Matplotlibでグラフ表示1(Matplotlibの基本確認) の章追加
    • Matplotlibでグラフ表示2(上下に2つのグラフ表示)
    • Matplotlibでグラフ表示3(パイチャート)
  • 2021/09/09 以下の章を追加
    • 新規PDFデータ追加
    • 主要なデータの確認及びグラフ表示

著者: NM Max

Created: 2021-09-09 木 22:55

Validate

自由研究2(厚労省のワクチンの成分がどこに行くのかを調べたPDF資料からの取り込みと、ソート)

自由研究2(厚労省のワクチンの成分がどこに行くのかを調べたPDF資料からの取り込みと、ソート)

1 概要

  • ファイザーワクチンの動物実験のデータがあり、どの臓器にどのくらい行くのかがわかる表がある。
    • PDFファイルで、処理しにくいのだが、多く行く臓器と、行く量がピークとなる時間の表を作成されている方がいたので、それをプログラムでやってみることにしました。

2 関連文書

2.1 厚労省の資料

2.2 使ったPythonパッケージ

3 厚労省のPDFファイルから2つの表を抜き出し、外部ファイルに保存

3.1 説明動画


3.2 厚労省のPDFファイルから2つの表を抜き出し、外部ファイルに保存する手順

3.2.1 仮想環境を用意し作業ディレクトリに入る

virtualenv covid-19
cd covid-19
source bin/activate

3.2.2 必要なパッケージインストール

pip install tabula-py
pip install ipython

3.2.3 元のPDFファイルをダウンロードして作業ディレクトリに保存

wget https://www.pmda.go.jp/drugs/2021/P20210212001/672212000_30300AMX00231_I100_1.pdf

3.2.4 PDFファイルのP.16,17から目的の表の一部を取り出して、ファイルに保存

  1. Pythonの対話環境起動
    ipython
    
  2. プログラムを動かして、表を作成して外部ファイルに保存
    import pandas as pd
    import tabula
    
    # pdfからデータ読み込み
    d0 = tabula.read_pdf("672212000_30300AMX00231_I100_1.pdf",pages="16,17")
    #d0 = tabula.read_pdf("672212000_30300AMX00231_I100_1.pdf",pages="16,17",stream=True)
    #d0 = tabula.read_pdf("672212000_30300AMX00231_I100_1.pdf",pages="16,17",lattice=True)
    
    # カラム名の変更
    d0[0].columns=['Sample','v1','v20','v21']
    d0[1].columns=['Sample','v1','v20','v21']
    
    for i,j in enumerate(d0[0]["v1"].isnull()):
      if j:
        d0[0]["Sample"][i-1]=d0[0]["Sample"][i-1]+" "+d0[0]["Sample"][i]
    for i,j in enumerate(d0[1]["v1"].isnull()):
      if j:
        d0[1]["Sample"][i-1]=d0[1]["Sample"][i-1]+" "+d0[1]["Sample"][i]
    
    
    # NaNのデータ削除
    d1=d0[0][d0[0]["v1"].notnull()]
    d2=d0[1][d0[1]["v1"].notnull()]
    
    # 経過時間名の配列作成
    idx=d2["v1"][1].replace(" h","h").split(" ")
    
    # 表の左側のみ残す
    d11=d1.drop([0,1]).drop(["v20","v21"],axis=1)
    d22=d2.drop([0,1]).drop(["v20","v21"],axis=1)
    
    # 列番号のふりなおし
    d111=d11.reset_index(drop=True)
    d222=d22.reset_index(drop=True)
    
    # データを配列に
    d1111=d111["v1"].str.rsplit(" ",expand=True)
    d2222=d222["v1"].str.rsplit(" ",expand=True)
    
    # 配列データが文字列なので、数値に変換
    for i in range(7):
      d1111[i] = pd.to_numeric(d1111[i])
      d2222[i] = pd.to_numeric(d2222[i])
    
    # 列の名前を経過時間に変更
    d1111.columns=idx
    d2222.columns=idx
    
    # 行毎の最高値と、そのインデックスを調べ、maxとtmaxという列を追加して格納
    tmpmax=d1111.max(axis=1)
    tmptmax=d1111.idxmax(axis=1)
    d1111["max"]=tmpmax
    d1111["tmax"]=tmptmax
    
    # 行毎の最高値と、そのインデックスを調べ、maxとtmaxという列を追加して格納2
    tmpmax=d2222.max(axis=1)
    tmptmax=d2222.idxmax(axis=1)
    d2222["max"]=tmpmax
    d2222["tmax"]=tmptmax
    
    # 一番左にenという列を追加して、英語名を入れる
    d1111.insert(0,"en",d111["Sample"])
    d2222.insert(0,"en",d222["Sample"])
    
    # P.16の表とp.17の表を結合してddという表にまとめる
    dd=pd.merge(d1111,d2222,how="outer")
    
    # 英語名を日本語に翻訳したk、jaという配列作成
    ja = ['脂肪組織', '副腎', '膀胱', '骨(大腿骨)', '骨髄(大腿骨)', '脳', '目', '心臓', '注射部位', '腎臓', '大腸', '肝臓', '肺', 'リンパ節(下顎)', 'リンパ節(腸間膜)', '筋', '卵巣(雌)', '膵臓', '脳下垂体', '前立腺(男性)', '唾液腺', '肌', '小腸', '脊髄', '脾臓', 'お腹', 'テスト(悪)', '胸腺', '甲状腺', '子宮(女性)', '全血', '血漿', '血液:血漿比']
    
    # ddの一番右に日本語名を"ja"という列名で追加
    dd.insert(len(dd.columns),"ja",ja)
    
    # pickleパッケージで、作成した表を外部ファイル seibun_data2.pickle に保存
    import pickle
    with open('seibun_data2.pickle', 'wb') as f:
        pickle.dump(dd,f)
    

3.2.5 作成した表データを表示してみる

                         en    0.25h       1h       2h       4h       8h      24h      48h      max   tmax         ja
0            Adipose tissue    0.057    0.100    0.126    0.128    0.093    0.084    0.181    0.181    48h       脂肪組織
1            Adrenal glands    0.271    1.480    2.720    2.890    6.800   13.800   18.200   18.200    48h         副腎
2                   Bladder    0.041    0.130    0.146    0.167    0.148    0.247    0.365    0.365    48h         膀胱
3              Bone (femur)    0.091    0.195    0.266    0.276    0.340    0.342    0.687    0.687    48h     骨(大腿骨)
4       Bone marrow (femur)    0.479    0.960    1.240    1.240    1.840    2.490    3.770    3.770    48h    骨髄(大腿骨)
5                     Brain    0.045    0.100    0.138    0.115    0.073    0.069    0.068    0.138     2h          脳
6                      Eyes    0.010    0.035    0.052    0.067    0.059    0.091    0.112    0.112    48h          目
7                     Heart    0.282    1.030    1.400    0.987    0.790    0.451    0.546    1.400     2h         心臓
8            Injection site  128.000  394.000  311.000  338.000  213.000  195.000  165.000  394.000     1h       注射部位
9                   Kidneys    0.391    1.160    2.050    0.924    0.590    0.426    0.425    2.050     2h         腎臓
10          Large intestine    0.013    0.048    0.093    0.287    0.649    1.100    1.340    1.340    48h         大腸
11                    Liver    0.737    4.630   11.000   16.500   26.500   19.200   24.300   26.500     8h         肝臓
12                     Lung    0.492    1.210    1.830    1.500    1.150    1.040    1.090    1.830     2h          肺
13  Lymph node (mandibular)    0.064    0.189    0.290    0.408    0.534    0.554    0.727    0.727    48h   リンパ節(下顎)
14  Lymph node (mesenteric)    0.050    0.146    0.530    0.489    0.689    0.985    1.370    1.370    48h  リンパ節(腸間膜)
15                   Muscle    0.021    0.061    0.084    0.103    0.096    0.095    0.192    0.192    48h          筋
16        Ovaries (females)    0.104    1.340    1.640    2.340    3.090    5.240   12.300   12.300    48h      卵巣(雌)
17                 Pancreas    0.081    0.207    0.414    0.380    0.294    0.358    0.599    0.599    48h         膵臓
18          Pituitary gland    0.339    0.645    0.868    0.854    0.405    0.478    0.694    0.868     2h       脳下垂体
19         Prostate (males)    0.061    0.091    0.128    0.157    0.150    0.183    0.170    0.183    24h    前立腺(男性)
20          Salivary glands    0.084    0.193    0.255    0.220    0.135    0.170    0.264    0.264    48h        唾液腺
21                     Skin    0.013    0.208    0.159    0.145    0.119    0.157    0.253    0.253    48h          肌
22          Small intestine    0.030    0.221    0.476    0.879    1.280    1.300    1.470    1.470    48h         小腸
23              Spinal cord    0.043    0.097    0.169    0.250    0.106    0.085    0.112    0.250     4h         脊髄
24                   Spleen    0.334    2.470    7.730   10.300   22.100   20.100   23.400   23.400    48h         脾臓
25                  Stomach    0.017    0.065    0.115    0.144    0.268    0.152    0.215    0.268     8h         お腹
26           Testes (males)    0.031    0.042    0.079    0.129    0.146    0.304    0.320    0.320    48h     テスト(悪)
27                   Thymus    0.088    0.243    0.340    0.335    0.196    0.207    0.331    0.340     2h         胸腺
28                  Thyroid    0.155    0.536    0.842    0.851    0.544    0.578    1.000    1.000    48h        甲状腺
29         Uterus (females)    0.043    0.203    0.305    0.140    0.287    0.289    0.456    0.456    48h     子宮(女性)
30              Whole blood    1.970    4.370    5.400    3.050    1.310    0.909    0.420    5.400     2h         全血
31                   Plasma    3.970    8.130    8.900    6.500    2.360    1.780    0.805    8.900     2h         血漿
32      Blood:Plasma ratioa    0.815    0.515    0.550    0.510    0.555    0.530    0.540    0.815  0.25h     血液:血漿比

3.2.6 表を読み込み、その臓器に行った量の多い順にソート、そして、結果をCSVで出力

import pickle
import pandas as pd
import tabula

# データ読み込み
with open('seibun_data2.pickle', 'rb') as f:
    dd=pickle.load(f)

print(dd.sort_values(["max","tmax"], ascending=False).reset_index(drop=True))
print()
print(dd.sort_values(["max","tmax"], ascending=False).reset_index(drop=True)[["en","max","tmax","ja"]])
print()
print(dd.sort_values(["max","tmax"], ascending=False).reset_index(drop=True)[["max","tmax","ja"]])

# csvファイルとして外部出力
dd.to_csv("seibun_out01.csv")
dd.to_csv("seibun_out02.csv", index=None)
  • 最後に出力した表
        max   tmax         ja
0   394.000     1h       注射部位
1    26.500     8h         肝臓
2    23.400    48h         脾臓
3    18.200    48h         副腎
4    12.300    48h      卵巣(雌)
5     8.900     2h         血漿
6     5.400     2h         全血
7     3.770    48h    骨髄(大腿骨)
8     2.050     2h         腎臓
9     1.830     2h          肺
10    1.470    48h         小腸
11    1.400     2h         心臓
12    1.370    48h  リンパ節(腸間膜)
13    1.340    48h         大腸
14    1.000    48h        甲状腺
15    0.868     2h       脳下垂体
16    0.815  0.25h     血液:血漿比
17    0.727    48h   リンパ節(下顎)
18    0.687    48h     骨(大腿骨)
19    0.599    48h         膵臓
20    0.456    48h     子宮(女性)
21    0.365    48h         膀胱
22    0.340     2h         胸腺
23    0.320    48h     テスト(悪)
24    0.268     8h         お腹
25    0.264    48h        唾液腺
26    0.253    48h          肌
27    0.250     4h         脊髄
28    0.192    48h          筋
29    0.183    24h    前立腺(男性)
30    0.181    48h       脂肪組織
31    0.138     2h          脳
32    0.112    48h          目

3.3 この章のまとめ

4 今後

  • 今後文書追加するかも、しないかも

5 この文書のチェンジログ

  • 2021/09/02 初版

著者: NM Max

Created: 2021-09-07 火 12:12

Validate

自由研究(新型コロナの遺伝子配列、PCR、ワクチンの遺伝子配列)





自由研究(新型コロナの遺伝子配列、PCR、ワクチンの遺伝子配列)

目次

1 概要

  • もともとは、PCR検査でワクチンの遺伝子は陽性になるのかを調べてみたくなって、やりはじめました
    • 最初PCRで調べてる遺伝子の場所を知らなかったので、それを確認するために以下を行いました。(PCR検査では新型コロナの1/300程度しか調べてません)
      • 最初中国の論文で新型コロナの遺伝子配列が調べられました。これを文字列化
        • ちなみに、1ヶ月かそれ以上前に変異種1万種類以上確認されてます。報道は少ない種類のように報道してますが2週間に1度くらいのペースで変異しまくり。
      • PCR検査のマニュアルが感染研にあるので、それに記載されている遺伝子配列が、中国で初めて確認された新型コロナの遺伝子配列のどこなのかを調べた
      • ファイザーワクチンの遺伝子配列が厚労省の文書に書かれていたので、それを文字列化(画像だったので、OCRで処理。正しく変換できてるかは面倒なので、未確認)
      • PCR検査の検査部位の遺伝子配列が、ファイザーワクチンの遺伝子配列に含まれているか確認した(結果含まれてないようです。)
        • 多少ちがってても増殖できることもあるらしいので、PCR検査で陽性になるかどうかを調べれたわけではありませんが、PCR検査で検査してる部分の遺伝子配列に完全に一致する遺伝子配列は存在してないと思われます。

2 関連文書

3 新型コロナ(一番最初の中国論文の)の遺伝子配列を文字列化して変数に保存

3.1 説明動画


3.2 新型コロナ(一番最初の中国論文の)の遺伝子配列を文字列化して変数に保存をおこなう手順

3.2.1 まず遺伝子配列が記載されているページをChomeで開き、外部ファイルに保存します。

3.2.2 BeautifulSoupが必要なので入ってなかったら入れる

  • Pipでインストールするなら多分以下
pip install beautifulsoup4
  • あるいはUbuntuでパッケージでインストールするなら以下のコマンド
sudo apt install python3-bs4

3.2.3 IPythonの起動

ipython3

or

ipython

3.2.4 データ処理

  • BeautifulSoupと正規表現用のreパッケージ読み込み
from bs4 import BeautifulSoup
import re
  • 先程保存したHTMLファイルを読み込む
soup = BeautifulSoup(open('Severe acute respiratory syndrome coronavirus 2 isolate Wuhan-Hu-1, co - Nucleotide - NCBI.html'), 'html.parser')
  • 遺伝子配列のある部分をクラス指定で取り出す
    • このクラス名で取り出せる事は元のHTMLをChromeで開き、デベロッパーツールで調べておいた
soup.find_all(text=True,class_="ff_line")
  • 上手く取り出せているので、1つめを文字列化してみる
soup.find_all(text=True,class_="ff_line")[0].text
  • 色々不必要な文字列が含まれてるので、小文字のアルファベット以外を除いてみる
re.sub(r"[^a-z]","",soup.find_all(text=True,class_="ff_line")[0].text)
  • mapを利用して、全部の要素に今の処理を実行
list(map(lambda x: re.sub(r"[^a-z]","",x.text),soup.find_all(text=True,class_="ff_line")))
  • 全部の要素に今の処理を行い結合したものを変数xxに保存
xx=''.join(list(map(lambda x: re.sub(r"[^a-z]","",x.text),soup.find_all(text=True,class_="ff_line"))))
  • 作成した遺伝子配列を表示してみる
print(xx)
  • 作成した遺伝子配列の長さを表示してみる
print("len(xx)")
print(len(xx))

3.3 この章のまとめ

  • 新型コロナウィルスの遺伝子配列(一番最初の中国論文)を文字列にして、変数に入れる処理を行なった

4 PCR検査マニュアルにある遺伝子配列を確認その1(感染研・地衛研専用のマニュアル)

4.1 説明動画


4.2 PCR検査マニュアルにある遺伝子配列を確認1 手順

4.2.1 IPythonを起動

ipython3

4.2.2 以下のコマンドを実行し、新型コロナウィルスの中国論文の遺伝子配列を変数に代入

from bs4 import BeautifulSoup
import re
soup = BeautifulSoup(open('Severe acute respiratory syndrome coronavirus 2 isolate Wuhan-Hu-1, co - Nucleotide - NCBI.html'), 'html.parser')
xx=''.join(list(map(lambda x: re.sub(r"[^a-z]","",x.text),soup.find_all(text=True,class_="ff_line"))))
print(len(xx))

4.2.3 逆順にする関数myrevと、ペアの遺伝子配列にする関数myconvを定義

def myrev(s):
    return ''.join(list(reversed(s)))

def myconv(s):
    return s.upper().replace("C","g").replace("G","c").replace("A","t").replace("T","a")

4.2.4 PCRの検査マニュアルに記載されている位置の遺伝子配列の確認(リアルタイムPCRの検査部分を確認) その1

  1. NIID-N2 セット 29125-29144 の範囲表示
    print(xx[29125-1:29144])
    
  2. NIID-N2 セット 29282-29263 の範囲表示
    print(xx[29263-1:29282])
    
    • これを逆順にしてペアの配列を表示
    print(myrev(myconv(xx[29263-1:29282])))
    
  3. NIID-N2 セット 29222-29241 の範囲表示
    print(xx[29222-1:29241])
    
  4. NIID-S ver.2 (S2)セット 24722-24742 の範囲表示
    print(xx[24722-1:24742])
    
  5. NIID-S ver.2 (S2)セット 24870-24851 の範囲表示
    print(xx[24851-1:24870])
    
    • これを逆順にしてペアの配列を表示
    print(myrev(myconv(xx[24851-1:24870])))
    
  6. NIID-S ver.2 (S2)セット 24793-24816 の範囲表示
    print(xx[24793-1:24816])
    

4.2.5 PCRの検査マニュアルに記載されている位置の遺伝子配列の確認(リアルタイムPCRの検査部分を確認) その2

  • re.finditerを使って探索
  1. NIID-N2 セット
    print(list(re.finditer("AAATTTTGGGGACCAGGAAC".lower(),xx)))
    print(list(re.finditer(myconv(myrev("TGGCACCTGTGTAGGTCAAC")).lower(),xx)))
    print(list(re.finditer("ATGTCGCGCATTGGCATGGA".lower(),xx)))
    
  2. NIID-S ver.2 (S2)セット
    print(list(re.finditer("CAGTCAGCACCTCATGGTGTA".lower(),xx)))
    print(list(re.finditer(myconv(myrev("AACCAGTGTGTGCCATTTGA")).lower(),xx)))
    print(list(re.finditer("TGCTCCTGCCATTTGTCATGATGG".lower(),xx)))
    

4.2.6 PCRの検査でチェックしている遺伝子配列の長さ

  • それぞれの長さの表示
print(list(map(lambda x:len(x),["AAATTTTGGGGACCAGGAAC","TGGCACCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA"])))
print(sum(list(map(lambda x:len(x),["AAATTTTGGGGACCAGGAAC","TGGCACCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA"]))))
print(list(map(lambda x:len(x),["CAGTCAGCACCTCATGGTGTA","AACCAGTGTGTGCCATTTGA","TGCTCCTGCCATTTGTCATGATGG"])))
print(sum(list(map(lambda x:len(x),["CAGTCAGCACCTCATGGTGTA","AACCAGTGTGTGCCATTTGA","TGCTCCTGCCATTTGTCATGATGG"]))))

4.3 この章のまとめ

  • 前の章の遺伝子配列で、PCR検査で検査している遺伝子配列を調べてみた(製品によって異る部分を調べているものがあります。)
  • 全長3万弱の遺伝子のほんのちょっとの部分しか調べていないことを確認

5 PCR検査マニュアルにある遺伝子配列を確認その2

5.1 説明動画


5.2 PCR検査マニュアルにある遺伝子配列を確認2 手順

5.2.1 ベースとなる新型コロナの遺伝子配列が旧バージョンなのでそれをダウンロード

5.2.2 IPythonを起動

ipython3

5.2.3 以下のコマンドを実行し、新型コロナウィルスの中国論文の遺伝子配列を変数に代入

  • xxにMN908947.3(最新版)
  • xx2にMN908947.1(旧バージョン)
from bs4 import BeautifulSoup
import re

def myrev(s):
    return ''.join(list(reversed(s)))

def myconv(s):
    return s.upper().replace("C","g").replace("G","c").replace("A","t").replace("T","a")

soup = BeautifulSoup(open('Severe acute respiratory syndrome coronavirus 2 isolate Wuhan-Hu-1, co - Nucleotide - NCBI.html'), 'html.parser')
xx=''.join(list(map(lambda x: re.sub(r"[^a-z]","",x.text),soup.find_all(text=True,class_="ff_line"))))
print(len(xx))

soup2 = BeautifulSoup(open('Wuhan seafood market pneumonia virus isolate Wuhan-Hu-1, complete geno - Nucleotide - NCBI.html'), 'html.parser')
xx2=''.join(list(map(lambda x: re.sub(r"[^a-z]","",x.text),soup2.find_all(text=True,class_="ff_line"))))
print(xx2)
print("len(xx2)")
print(len(xx2))

5.2.4 探索部分を関数に定義

def mycheck(xx,x1="CACATTGGCACCCGCAATC",x2="GAGGAACGAGAAGAGGCTTG",x3="ACTTCCTCAAGGAACAACATTGCCA"):
  s=xx.upper()
  return(list(re.finditer(x1.upper(),s)),list(re.finditer(myconv(myrev(x2.upper())).upper(),s)) ,list(re.finditer(x3.upper(),s)))

def mycheck2(xx):
  return mycheck(xx,"AAATTTTGGGGACCAGGAAC","TGGCAGCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA")

5.2.5 それぞれの長さの表示

  • Nセット
print(list(map(lambda x:len(x),["CACATTGGCACCCGCAATC","GAGGAACGAGAAGAGGCTTG","ACTTCCTCAAGGAACAACATTGCCA"])))
print(sum(list(map(lambda x:len(x),["CACATTGGCACCCGCAATC","GAGGAACGAGAAGAGGCTTG","ACTTCCTCAAGGAACAACATTGCCA"]))))
  • N2セット
print(list(map(lambda x:len(x),["AAATTTTGGGGACCAGGAAC","TGGCAGCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA"])))
print(sum(list(map(lambda x:len(x),["AAATTTTGGGGACCAGGAAC","TGGCAGCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA"]))))

5.2.6 それぞれの長さの全体にしめる割合を確認

  • Nセット
# Nセット
print(len(xx2)/sum(list(map(lambda x:len(x),["CACATTGGCACCCGCAATC","GAGGAACGAGAAGAGGCTTG","ACTTCCTCAAGGAACAACATTGCCA"]))))
  • N2セット
# N2セット
print(len(xx2)/sum(list(map(lambda x:len(x),["AAATTTTGGGGACCAGGAAC","TGGCAGCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA"]))))
  • NセットとN2セット
# NセットとN2セット
print(len(xx2)/(sum(list(map(lambda x:len(x),["AAATTTTGGGGACCAGGAAC","TGGCAGCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA"])))+sum(list(map(lambda x:len(x),["AAATTTTGGGGACCAGGAAC","TGGCAGCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA"])))))

5.2.7 PCRの検査マニュアルに記載されている位置の遺伝子配列の確認(リアルタイムPCRの検査部分を確認) その1

  • N1セットのチェック
    • なぜか、1つ場所がズレてる。原因は不明。遺伝子配列は存在してる
# 28723-28741
# 28850-28831
# 28770-28794
mycheck(xx2)
  • N2セットのチェック
    • なぜか、1つ場所がズレてる。原因は不明。遺伝子配列は存在してる
# 29142-29161
# 29299-29280
# 29239-29258
mycheck2(xx2)

5.3 この章のまとめ

  • 通常のPCR検査のマニュアルにあるチェック部分の遺伝子配列を調べた

6 ファイザーワクチンのmRNAの遺伝子配列を文字列化して変数に格納

6.1 説明動画


6.2 ファイザーワクチンのmRNAの遺伝子配列を文字列化して変数に格納 手順

6.2.1 ファイザーの遺伝子配列情報が記載されてる書類をダウンロード

  • 審議結果報告書(ファイザーの遺伝子配列情報が記載されてる) https://www.mhlw.go.jp/content/10601000/000739089.pdf
  • Yが特殊で”一メチルシュウドウウリジン”というのがいれられていて、通常のmRNAではない。長寿命になるように工夫されている。(報道のように普通のmRNAじゃなく、長寿命化)

6.2.2 PDFファイルから、遺伝子情報のあるページを抜き出し

  • 000739089-04.ppm 000739089-05.ppm 000739089-06.ppm のファイルが生成される。画像フォーマットはppm
pdftoppm -l 6 -f 4 000739089.pdf 000739089

6.2.3 画像ファイルから余計な部分を除去

  • 私はGIMPを用いて切り抜きしました。他の画像処理ツールでもOK

6.2.4 OCRソフト( tesseract )で文字列化

  • psmオプションつけないと、列毎の文字列になってしまい、後処理が面倒になる
  • 000739089-04.txt 000739089-05.txt 000739089-06.txt という3つのテキストファイルが生成される。
tesseract -l eng  --psm 4 000739089-04.png 000739089-04
tesseract -l eng  --psm 4 000739089-05.png 000739089-05
tesseract -l eng  --psm 4 000739089-06.png 000739089-06

6.2.5 Python起動

ipython3

or

ipython

6.2.6 OCRソフト( tesseract )で文字列化

import re
yy=""
with open("000739089-04.txt") as f:
    yy+=re.sub(r"[^a-zA-Z]","",f.read())
with open("000739089-05.txt") as f:
    yy+=re.sub(r"[^a-zA-Z]","",f.read())
with open("000739089-06.txt") as f:
    yy+=re.sub(r"[^a-zA-Z]","",f.read())

6.2.7 文字列の確認

  1. 文字列を表示してみる(確認)
    print(yy)
    
  2. 文字列に含まれている文字を確認(確認)
    print(set(yy))
    
  3. 文字列の長さを確認(確認)
    • 文書によると4284
    # 4284
    print(len(yy))
    

6.3 この章のまとめ

  • ファイザーワクチンのmRNAの遺伝子配列を文字列化

7 ファイザーワクチンのmRNAの遺伝子配列にPCR検査している遺伝子配列が含まれているか確認

7.1 説明動画


7.2 ファイザーワクチンのmRNAの遺伝子配列にPCR検査している遺伝子配列が含まれているか確認 手順

7.2.1 Pythonの対話環境起動

ipython3

or

ipython

7.2.2 新型コロナの遺伝子配列の読み込み,今まで使った関数の定義

from bs4 import BeautifulSoup
import re

soup = BeautifulSoup(open('Severe acute respiratory syndrome coronavirus 2 isolate Wuhan-Hu-1, co - Nucleotide - NCBI.html'), 'html.parser')
xx=''.join(list(map(lambda x: re.sub(r"[^a-z]","",x.text),soup.find_all(text=True,class_="ff_line"))))
print(len(xx))

soup2 = BeautifulSoup(open('Wuhan seafood market pneumonia virus isolate Wuhan-Hu-1, complete geno - Nucleotide - NCBI.html'), 'html.parser')
xx2=''.join(list(map(lambda x: re.sub(r"[^a-z]","",x.text),soup2.find_all(text=True,class_="ff_line"))))

yy=""
with open("000739089-04.txt") as f:
    yy+=re.sub(r"[^a-zA-Z]","",f.read())
with open("000739089-05.txt") as f:
    yy+=re.sub(r"[^a-zA-Z]","",f.read())
with open("000739089-06.txt") as f:
    yy+=re.sub(r"[^a-zA-Z]","",f.read())

def myrev(s):
    return ''.join(list(reversed(s)))

def myconv(s):
    return s.upper().replace("C","g").replace("G","c").replace("A","t").replace("T","a")

def mycheck(xx,x1="CACATTGGCACCCGCAATC",x2="GAGGAACGAGAAGAGGCTTG",x3="ACTTCCTCAAGGAACAACATTGCCA"):
  s=xx.upper()
  return(list(re.finditer(x1.upper(),s)),list(re.finditer(myconv(myrev(x2.upper())).upper(),s)) ,list(re.finditer(x3.upper(),s)))

def mycheck1(xx):
  return mycheck(xx,x1="CACATTGGCACCCGCAATC",x2="GAGGAACGAGAAGAGGCTTG",x3="ACTTCCTCAAGGAACAACATTGCCA")

def mycheck2(xx):
  return mycheck(xx,"AAATTTTGGGGACCAGGAAC","TGGCAGCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA")

def mycheck01(xx):
  return mycheck(xx,"AAATTTTGGGGACCAGGAAC","TGGCACCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA")

def mycheck02(xx):
  return mycheck(xx,"CAGTCAGCACCTCATGGTGTA","AACCAGTGTGTGCCATTTGA","TGCTCCTGCCATTTGTCATGATGG")
  • check1(通常のPCR検査)とcheck02(研究所用)は同じ遺伝子配列を調べているようです。

7.2.3 チェック関数の動作確認

print(mycheck01(xx))
print(mycheck02(xx))
print(mycheck1(xx))
print(mycheck2(xx))

7.2.4 ワクチンの遺伝子配列のYの文字をTに変換

  • Yが特殊で”一メチルシュウドウウリジン”なので、Tに変更して別変数に入れる
yyy=yy.replace("Y","T")

7.2.5 ワクチンの遺伝子配列にPCR検査で調べている遺伝子配列があるかチェック

print(mycheck01(yyy))
print(mycheck02(yyy))
print(mycheck1(yyy))
print(mycheck2(yyy))
  • 逆順にして確認
print(mycheck01(myrev(yyy)))
print(mycheck02(myrev(yyy)))
print(mycheck1(myrev(yyy)))
print(mycheck2(myrev(yyy)))

7.3 チェック結果

  • ファイザーワクチンのmRNAの遺伝子配列内にPCR検査で調べている遺伝子配列は見つからなかった
    • 多少ずれてても増幅すると聞いているので、このチェックがPCR検査で陽性にならないかどかのチェックにはなってないと思います。
    • 新型コロナの遺伝子の一部が含まれてると聞いていたのに、なぜ含まれないのか、理由は不明。
    • PCR検査で陽性にならないように、同義の遺伝子配列にする等の工夫がなされてるのかも?そこ調べてません。
      • 同義置換して一致するか確認してみるのも面白そう
    • ワクチンの遺伝子配列が正しくOCRで変換できてるかは、面倒だったので未確認なので、そこらあるかも

7.4 この章のまとめ

  • ファイザーワクチンのmRNAの遺伝子配列内にPCR検査で調べている遺伝子配列そのものは見つからなかった

8 PCR検査で調べている遺伝子配列が 新型コロナの MN908947.3(最新バージョン),MN908947.1(昔のバージョン)に存在するか確認

8.1 説明動画


8.2 PCR検査で調べている遺伝子配列が 新型コロナの MN908947.3(最新バージョン),MN908947.1(昔のバージョン)に存在するか確認 の手順

  • 前の章のプログラムを動かして準備完了

8.2.1 チェック関数の動作確認

  1. xxにMN908947.3(最新版)も前の章やったが、再度やっておく
    print(mycheck01(xx))
    print(mycheck02(xx))
    print(mycheck1(xx))
    print(mycheck2(xx))
    
  2. xx2にMN908947.1(旧バージョン)
    • check01の2番目のTGGCACCTGTGTAGGTCAACがマッチしてない。MN908947.1(旧バージョン)にはこの遺伝子配列そのものは無いようだ。
    print(mycheck01(xx2))
    print(mycheck02(xx2))
    print(mycheck1(xx2))
    print(mycheck2(xx2))
    

8.3 この章のまとめ

  • 研究所でのPCR検査で用いる遺伝子配列群1の終了用の遺伝子配列が、新型コロナの遺伝子配列 MN908947.1(旧バージョン)に含まれていない事を確認した。
    • 研究所用のマニュアルの遺伝子配列をこれに選択した理由があるはずだが、理由は不明。
  • 通常検査用のプライマー配列そのものは、新型コロナの遺伝子配列 MN908947.3(最新バージョン)には含まれていない。しかし、マニュアルには検出感度に問題無いとの記載あり。
    • 遺伝子配列が多少異なっても増幅するとは聞いていたが、こちらも最新バージョンの遺伝子配列をベースにせずに、旧バージョンの遺伝子配列からプライマーを決めた理由は良くわからない。

9 ファイザーワクチンのmRNAの遺伝子配列にPCR検査している遺伝子配列に似てる遺伝子配列が含まれているか確認1

  • 審議結果報告書(ファイザーの遺伝子配列情報が記載されてる) https://www.mhlw.go.jp/content/10601000/000739089.pdf
  • PCR検査で調べている遺伝子配列と1文字違いの遺伝子配列がファイザーの遺伝子配列に含まれているか調べてみる

9.1 説明動画


9.2 ファイザーワクチンのmRNAの遺伝子配列にPCR検査している遺伝子配列に似てる遺伝子配列が含まれているか確認1 の手順

  • Pythonの対話環境を起動して以下を実行
from bs4 import BeautifulSoup
import re

# 新型コロナの遺伝子配列読み込み 変数xx MN908947.3 最新バージョン
soup = BeautifulSoup(open('Severe acute respiratory syndrome coronavirus 2 isolate Wuhan-Hu-1, co - Nucleotide - NCBI.html'), 'html.parser')
xx=''.join(list(map(lambda x: re.sub(r"[^a-z]","",x.text),soup.find_all(text=True,class_="ff_line"))))
print(len(xx))

# 新型コロナの遺伝子配列読み込み 変数xx2 MN908947.1 旧タイプ
soup2 = BeautifulSoup(open('Wuhan seafood market pneumonia virus isolate Wuhan-Hu-1, complete geno - Nucleotide - NCBI.html'), 'html.parser')
xx2=''.join(list(map(lambda x: re.sub(r"[^a-z]","",x.text),soup2.find_all(text=True,class_="ff_line"))))

# ファイザーワクチンの遺伝子配列を 変数yyに
yy=""
with open("000739089-04.txt") as f:
    yy+=re.sub(r"[^a-zA-Z]","",f.read())
with open("000739089-05.txt") as f:
    yy+=re.sub(r"[^a-zA-Z]","",f.read())
with open("000739089-06.txt") as f:
    yy+=re.sub(r"[^a-zA-Z]","",f.read())

# ファイザーワクチンの遺伝子配列を 変数yyyに("一メチルシュウドウウリジン" YをTに変換)
yyy=yy.replace("Y","T")

# 関数群
# 文字列を逆順に変換
def myrev(s):
    return ''.join(list(reversed(s)))

# CとG、AとTを入替え
def myconv(s):
    return s.upper().replace("C","g").replace("G","c").replace("A","t").replace("T","a")

# チェック関数
def mycheck(xx,x1="CACATTGGCACCCGCAATC",x2="GAGGAACGAGAAGAGGCTTG",x3="ACTTCCTCAAGGAACAACATTGCCA"):
  s=xx.upper()
  return(list(re.finditer(x1.upper(),s)),list(re.finditer(myconv(myrev(x2.upper())).upper(),s)) ,list(re.finditer(x3.upper(),s)))

# チェック関数1 通常のPCR検査で調べる遺伝子配列群1
def mycheck1(xx):
  return mycheck(xx,x1="CACATTGGCACCCGCAATC",x2="GAGGAACGAGAAGAGGCTTG",x3="ACTTCCTCAAGGAACAACATTGCCA")

# チェック関数2 通常のPCR検査で調べる遺伝子配列群2
def mycheck2(xx):
  return mycheck(xx,"AAATTTTGGGGACCAGGAAC","TGGCAGCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA")

# チェック関数01 研究所でのPCR検査で用いる遺伝子配列群1
def mycheck01(xx):
  return mycheck(xx,"AAATTTTGGGGACCAGGAAC","TGGCACCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA")

# チェック関数02 研究所でのPCR検査で用いる遺伝子配列群2
def mycheck02(xx):
  return mycheck(xx,"CAGTCAGCACCTCATGGTGTA","AACCAGTGTGTGCCATTTGA","TGCTCCTGCCATTTGTCATGATGG")

# 元の文字列 s の index 番目の文字列を newstringに置き換える関数 
def myreplace(s,index,newstring):
  return s[:index] + newstring + s[index + 1:]

# 元の文字列 s の1文字を"."に置き換えて出来る文字列を生成する関数
def mymake1replace(s):
  ans=[]
  for i in range(len(s)):      
      ans.append(myreplace(s,i,"."))
  return ans

# mycheckv2で利用する関数
def bmycheckv2(tt,ss):
  return (tt,list(re.finditer(tt,ss)))

# mycheckv2で利用する関数
def brmycheckv2(tt,ss):
  ttt=myrev(myconv(tt).upper())
  return (ttt,list(re.finditer(ttt,ss)))

# 遺伝子配列sに tという遺伝子配列、tの1文字違いの文字列が含まれているか検索する関数
def mycheckv2(s,t):
  ss=s.upper()
  tt=t.upper()
  ans1=[]
  ans2=[]
  for i in mymake1replace(tt):
    t1,t2 = bmycheckv2(i,ss)
    if len(t2)!=0: 
      ans1.append((t1,t2))
    tt1,tt2 = brmycheckv2(i,ss)
    if len(tt2)!=0: 
      ans1.append((tt1,tt2))
  return (ans1,ans2)

#ヒットしない文字列を filter


cs1=("CACATTGGCACCCGCAATC","GAGGAACGAGAAGAGGCTTG","ACTTCCTCAAGGAACAACATTGCCA")
cs2=("AAATTTTGGGGACCAGGAAC","TGGCAGCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA")
cs3=("AAATTTTGGGGACCAGGAAC","TGGCACCTGTGTAGGTCAAC","ATGTCGCGCATTGGCATGGA")
cs4=("CAGTCAGCACCTCATGGTGTA","AACCAGTGTGTGCCATTTGA","TGCTCCTGCCATTTGTCATGATGG")

### 動作確認
# 研究所用の検査マニュアル 新型コロナの遺伝子配列読み込み 変数xx MN908947.3 最新バージョン
for i in cs3:
  print("########## "+i)
  print(mycheckv2(xx,i))
  print()
  print()

for i in cs4:
  print("########## "+i)
  print(mycheckv2(xx,i))
  print()
  print()

# 通常の検査マニュアル 新型コロナの遺伝子配列読み込み 変数xx2 MN908947.1 旧バージョン
for i in cs1:
  print("########## "+i)
  print(mycheckv2(xx2,i))
  print()
  print()

for i in cs2:
  print("########## "+i)
  print(mycheckv2(xx2,i))
  print()
  print()

# 2番目のプライマーが発見出来ない組み合わせ # cs2 通常の検査用マニュアル
for i in cs2:
  print("########## "+i)
  print(mycheckv2(xx,i))
  print()
# 
# 2番目のプライマーが発見出来ない組み合わせ # cs3 研究所用のマニュアル

for i in cs3:
  print("########## "+i)
  print(mycheckv2(xx2,i))
  print()


for i in cs1:
    print("########## "+i)
    print(mycheckv2(yyy,i))
    print()
    print()

for i in cs2:
    print("########## "+i)
    print(mycheckv2(yyy,i))
    print()
    print()

for i in cs3:
    print("########## "+i)
    print(mycheckv2(yyy,i))
    print()
    print()

for i in cs4:
    print("########## "+i)
    print(mycheckv2(yyy,i))
    print()
    print()

9.3 この章のまとめ

  • PCR検査で調べている遺伝子配列と1文字違いの遺伝子配列もファイザーの遺伝子配列から見つけることは出来なかった。
    • 元のOCRでの変換が正しいとして

10 今後

  • 今後文書追加するかも、しないかも

11 この文書のチェンジログ

  • 2021/08/16 初版
  • 2021/08/17 PCR検査マニュアルにある遺伝子配列を確認その1(感染研・地衛研専用のマニュアル) の章,PCR検査マニュアルにある遺伝子配列を確認その2 の章を追加
  • 2021/08/18 PCR検査マニュアルにある遺伝子配列を確認その2 の章を修正
  • 2021/08/18 ファイザーワクチンのmRNAの遺伝子配列を文字列化して変数に格納 の章追加
  • 2021/08/18 ファイザーワクチンのmRNAの遺伝子配列にPCR検査している遺伝子配列が含まれているか確認 の章追加
  • 2021/09/03 PCR検査で調べている遺伝子配列が 新型コロナの MN908947.1(昔のバージョン)に存在するか確認 の章, ファイザーワクチンのmRNAの遺伝子配列にPCR検査している遺伝子配列に似てる遺伝子配列が含まれているか確認1 の章追加

著者: NM Max

Created: 2021-09-05 日 09:21

Validate