【備忘録】pythonのdatetimeモジュールのクラス種類と使い方

pythonのdatetimeモジュールで利用できるクラスの種類と違いを実例で確認する。

公式サイト

クラスの種類

6種類のクラス（データ型）が用意されている

1. datetime.date：日付(年、月、日)
2. datetime.time：時刻(時、分、秒、マイクロ秒)
3. datetime.datetime：日付＋時刻
4. datetime.timedelta：時間差分
5. datetime.timezone：タイムゾーンの設定
6. datetime.tzinfo：タイムゾーン情報

各クラスの中身

0. 補足条件

datetimeモジュールはdtとして呼び出す。

datetimeモジュールの呼び出し

import datetime as dt

▼実際の使用例
datetime.dateはdt.dateとして記述

※datetime.datetime.dateはエラーになる
（name 'datetime' is not defined）

1. datetime.date

dt.date.today()の実行結果
datetime.date(2020, 3, 18)
　└ datetime.date型
　 └ (2020, 3, 18)＝(年, 月, 日)

⇒ 年、月、日の3つのデータが取得できる。
(属性：year, month, day)

2. datetime.time

dt.timeでメソッドは使えなさそう。。
dt.dateのように、現在時刻を取得しようと思ってもできない。

• dt.time.now()： エラー
• dt.time.today() ：エラー

type object 'datetime.time' has no attribute 'today'

■使い方
時刻を自分で指定する。

日付データが不要で、時刻データのみ呼び出したり、変更したい場合に有効。

• datetime.time型がもつ、hours, minute, second, microsecondを個別に指定。
• 指定したデータを呼び出す。（例：.hour）
• 指定したデータを変更する （例：.replace(hour=1)）

指定した数値

時刻の指定と呼び出し（datetime.time）

#時刻を設定
timeA = dt.time(10, 30, 45, 123456) 

#時刻を呼び出す
timeA.hour  #実行結果：10
timeA.minute #実行結果：30
timeA.second #実行結果：45
timeA.microsecond #実行結果：123456

時刻の変更（replace()）

#時刻を設定
timeA = dt.time(10, 30, 45, 123456) 

#(1)個別に変更する
timeA = timeA.replace(hour=5)  
timeA

#実行結果
# →　datetime.time(5, 30, 45, 123456) 


#(2)まとめて変更する
timeA = timeA.replace(hour=1, minute=2, second=30, microsecond=400)
timeA

#出力結果
# → datetime.time(1, 2, 30, 400)

3. datetime.datetime

dt.datetime.today()の実行結果
datetime.datetime(2020, 3, 18, 7, 42, 54, 95450)
　└ datetime.datetime型
　 └ (2020, 3, 18, 7, 42, 54, 95450)
　　 (年, 月, 日,時,分,秒,マイクロ秒)

⇒ 年、月、日、時、分、秒、マイクロ秒の7つのデータが取得できる。
(属性：year, month, day, hour, minute, microsecond)

※裏側ではtzinfo（tz：タイムゾーン情報）も所持している。

▼tzinfo（tz：タイムゾーン情報）
カッコなしのnowなどでメソッドを実行した場合に情報が出力される。

タイムゾーン情報を表示する

dt.datetime.now
<function datetime.now(tz=None)>

■now()とtoday()
datetime.datetimeでよく使うメソッド。
※どちらも同じ

よく使うメソッドnowとtoday

#.now()
dtNow = dt.datetime.now()
dtNow  #datetime.datetime(2020, 3, 18, 21, 28, 13, 409431)

#.today()
dtToday = dt.datetime.today()
dtToday #datetime.datetime(2020, 3, 18, 22, 45, 10, 518281)

■combine(A, B)
日付の型と時刻の型を組み合わせることが可能。
datetime.datetime.combine(A, B)

よく使うメソッドcombine

#.combine(A, B)
dateA = dt.date(2020,5,6)
timeA = dt.time(10,20,0)

dt.datetime.combine(dateA, timeA)

#出力結果
# → datetime.datetime(2020, 5, 6, 10, 30, 0)

※combine(A, B)はA= datetime.date型、B= datetime.time型と決まっている。

cobmbine(timeA, dateB)など、指定の型と異なる場合はエラーになる。

combine() argument 1 must be datetime.date, not datetime.time

■replace()
引数で指定した数値を変更できる（複数指定可能）

replace

#.replace()

dtNow =dt.datetime.now() 
dtNow # datetime.datetime(2020, 3, 18, 22, 10, 15, 517216)

dtNow = dtNow.replace(year=2018)
dtNow = dtNow.replace(hour=5, minute=30, second=30)
dtNow

#出力結果
# → datetime.datetime(2018, 3, 18, 5, 30, 30, 474609)

■曜日を出力する
weekday()　：0～6（月曜日0～日曜日6）
isoweekday()　：1～6（月曜日1~日曜日7）

出力	月	火	水	木	金	土	日
weekday()	0	1	2	3	4	5	6
isoweekday()	1	2	3	4	5	6	7

曜日の出力

today = dt.datetime.today()  #datetime.datetime(2020, 3, 18, 22, 30, 46, 646567)

#.weekday()
today.weekday() #出力：2　→ 水曜日

#.isoweekday()
today.weekday() #出力：3　→ 水曜日

■日付 or 時刻のみ出力する
date()　：year、month、day
time()　：hour、minute、second、microsecond

日付・時刻の出力

now = dt.datetime.now()
#datetime.datetime(2020, 3, 18, 22, 35, 50, 177279)

#.time()
now.time() #出力：datetime.time(22, 38, 1, 997649)

#.date()
now.date() #datetime.date(2020, 3, 18)

4. datetime.timedelta

・日付や時間の足し算・引き算ができる。
・日付同士をひいた結果がtimedeltaで表示される。

• timedeltaは足し算・引き算処理をするものではい。
• ★どの項目をいくつ足すか・引くかの指定
  • dt.timedelta(days=2) ：2日足す
  • dttimedelta(hours=-5) :5時間引く
• 指定できる引数(7種類。※複数形
  • weeks （：7daysの倍数)
  • days
  • hours
  • minutes
  • seconds
  • milliseconds
  • microseconds
• 指定できない引数（年、月は不可）
  • years
  • months
• デフォルトの引数は0

指定できない引数や、間違った引数(単数形など)はエラーになる

'years' is an invalid keyword argument for __new__()

使い方
(1)2020年1月1日に2週間足す
+ dt.timedelta(weeks=2)

timedeltaで足し算

newyear = dt.datetime(2020, 1, 1)  #datetime型で2020年1月1日を定義
newyear + dt.timedelta(weeks=2)

#実行結果
#datetime.datetime(2020, 1, 15, 0, 0)
#　→ 2020年1月15日

(2)2020年1月1日から25時間引く
+ dt.timedelta(hours=-25)

timedeltaで引き算

newyear = dt.datetime(2020, 1, 1)  #datetime型で2020年1月1日を定義
newyear + dt.timedelta(hours=-25)

#実行結果
datetime.datetime(2019, 12, 30, 23, 0)
#　→ 2019年12月30日23時

補足
(3)所持していないデータ部分も計算される

▼実例
・date型(時間データを持たない)に時間足してもエラーにならない。
・結果は計算されたものが表示される
　　└ 非表示部で計算している

型が所持していないデータの取扱

newyear2 = dt.date(2020, 1, 1)　　#datetime.date型（時間データを持たない）
newyear2 + dt.timedelta(hours=100) #100時間足す

#実行結果
#datetime.date(2020, 1, 5)
# → 2020年1月5日

(3)日付の差分の結果
ある日付同士の間隔を算出する。

• 使える型はdatetime.datetimeとdatetime.date
• datetime.timeは使えない
• 比較する方は同じ必要がある
  • datetime.datetimeとdatetime.dateではエラーになる
• 足し算はできない（エラーになる）

2つの日付の間隔を求める（事例①）

dateA = dt.date(2020,1,1)
dateB = dt.date(2020,1,9)

dateA - dateB  #出力：datetime.timedelta(days=-8)
dateB - dateA  #出力：datetime.timedelta(days=8)

2つの日付の間隔を求める(事例②）

pastA = dt.datetime(2000, 1, 1)
now = dt.datetime.now()  #出力：datetime.datetime(2020, 3, 18, 22, 56, 50, 604259)

now - pastA

#出力結果
#datetime.timedelta(days=7382, seconds=82358, microseconds=594717)

※datetime.timedeltaにyears, monthsはない

足し算はエラー(※「+」はサポートしてない）

unsupported operand type(s) for +

5. datetime.timezone

自分でタイムゾーンを設定する場合に使用。

dt.timezone(dt.timedelta(hours=AAA), 'BBB')

設定する値は2つ。

• AAA：基準からどれだけズラすか
  • +/-0~24 ※24以下
    • 24はエラー
• BBB：設定したタイムゾーンにつける名称
  • 呼び出すときに引数で使う

▼設定方法
①イギリス(UCT。基準：+0時間)　[≒GMT]
datetime.timezone(datetime.timedelta(0), 'UCT')

②トロント（EDT。-4）
datetime.timezone(datetime.timedelta(-4), 'EST')

③日本（JST。+9）
dt.timezone(dt.timedelta(hours=+9), 'JST')

▼実際の使い方の例
トロントの時刻(EDT)を設定し、現在時刻を呼び出す。

timezoneの使い方(トロントの時刻)

EDT = dt.timezone(dt.timedelta(hours=-4),'EDT')
EDT　#出力：datetime.timezone(datetime.timedelta(days=-1, seconds=72000), 'EDT')

dt.datetime.now(EDT)

#出力結果
# datetime.datetime(2020, 3, 18, 10, 34, 8, 685415, tzinfo=datetime.timezone(datetime.timedelta(days=-1, seconds=72000), 'EDT'))

■補足:UTCとGMTの違い

①UTC：協定世界時間（Coordinated Universal Time)
　└ 現在の世界標準時間
　└ セシウム原子の振動数から算出

②GMT：グリニッジ標準時間（Greenwich Mean Time）
　└ 過去の標準時間
　└ グリニッジ天文台の計測結果から算出
　└ グリニッジはロンドンのすぐ南の都市

• どちらも大して変わらない（UTC≒GMT）
• イギリスの時間
  • イギリスが基準：hours=±０

6. datetime.tzinfo

タイムゾーンや夏時間など、カスタマイズ可能な時刻修正の概念を提供する（とのこと）。

はじめに

サブスクリプションモデル系はデータが豊富になりやすい（使い放題だから使う）と思っています。
私の場合は最近、音楽をものすごく聞くため、Spotify APIを使用すると面白いかなと思いました。

なので、今回はSpotify APIを使用して何かアプリケーションを作ってみます。

アプリケーションを登録しよう

https://developer.spotify.com/dashboard/

ここから登録可能です。
すでにアカウントがあればそれでログインすればすぐ使えるようになります。

登録中に商用利用するかというのも聞かれますので、正確に答えましょう。
アプリケーションを作成したら、EDIT SETTINGSからOAuthの情報を入力可能となります。
Webだったりスマホだったり色々と選択肢があります。
自分に合ったものを登録しましょう。

今回は後述しますが、厳密にWebサービスから叩くわけではないため、適当なURLを設定しています。

認証周りを作る

https://developer.spotify.com/documentation/general/guides/authorization-guide/

詳しくは上記に記載があります。
今回はWebアプリケーションとして作成しますが、実態としてはGoogle Cloud Functionsから実行するため、Refresh Tokenが必要となります。
Spotifyの仕様上、Access Tokenは1時間、Refresh Tokenは何かしら無効化するまで有効であり続けるとなっています。
なので、今回はアプリケーション作成前にRefresh Tokenを手で作って、コード上でAccess Tokenへ変換 -> SpotifyのAPIを叩くという手順を取ります。

1. code取得

まずはcodeの取得を行うためのURLを生成します。
詳細は省きますが、最終的に以下のようなURLになります。
client_idとredirect_uriは環境に合わせてください。

また、Spotify APIはscopeがとても細かく区切られています。
下記に記載されていますので、必要なscopeを指定するようにする必要があります。
https://developer.spotify.com/documentation/general/guides/scopes/

https://accounts.spotify.com/authorize?client_id=${CLIENT_ID}&response_type=code&redirect_uri=${REDIRECT_URL}&scope=user-read-private%20user-read-email&state=34fFs29kd09

出来上がったら適当なブラウザに貼り付けてEnter押すと設定したredirect_uriに帰ってきます。
ここでエラーになっていても特に問題はありません。
URL自体は変わっているため、そのURLからcodeの部分だけを次で使用します。

2. Refresh Tokenの取得

次はPOSTを投げないといけません。
なので、curlで投げます。
具体的には下記コマンドを取れたcodeと環境に合わせたclient_id、client_secretを設定して投げます。
問題なければトークンが返ってきますので、Refresh Tokenを保持しておきます。

curl --data "code=${CODE}" --data "client_id=${CLIENT_ID}" --data "client_secret=${CLIENT_SECRET}" --data "redirect_uri=http://localhost/callback" --data "grant_type=authorization_code" https://accounts.spotify.com/api/token

from Refresh Token to Access Token

Refresh TokenからAccess Tokenへの変換はAPI叩けば一発です。
その際はHeaderにAuthorization: Basic <base64 encoded client_id:client_secret>を含めないといけなく、それを作るのが面倒くさいと思ったぐらいであとは簡単です。

一応、自分が作ったpythonコードを貼っておきます。

from dotenv import load_dotenv

import os
import requests
import base64
import json

load_dotenv(verbose=True)

client_id = os.environ.get("SPOTIFY_CLIENT_ID", default="")
client_secret = os.environ.get("SPOTIFY_CLIENT_SECRET", default="")
token = base64.b64encode((client_id + ":" + client_secret).encode("utf-8")).decode(
    "utf-8"
)

headers = {"Authorization": "Basic " + token}
data = {
    "refresh_token": os.environ.get("REFRESH_TOKEN"),
    "grant_type": "refresh_token",
}
response = requests.post(
    "https://accounts.spotify.com/api/token", data=data, headers=headers
)  # noqa: E501
print(json.loads(response.text)["access_token"])

APIコールする部分を作成する

APIコールは取れたAccess TokenをBearerに設定して投げるだけです。
下に自分が作ろうとしているランキングデータ生成についてのサンプルコードを載せます。

    header = {"Authorization": "Bearer " + access_token}
    data = {
        "limit": 50,
        "time_range": "short_term",
    }
    response = requests.get(
        "https://api.spotify.com/v1/me/top/tracks", params=data, headers=header
    )

終わりに

ということでこれで無事Spotify API使ってアプリケーションが作れるようになりました。
Spotifyはどこに何があるかが分かりにくい反面、API仕様については細かく書かれています。
出来ることも意外と多いので、API作成の初歩としては中々良いサンプルとして使えるというのが感想です。

ちなみに今作ってるのはTwitterに月に1回、前月の聞いてたランキングを投稿するというものです。
まだ投稿部分が出来ていないため、これから作ろうと思います。

他のOAuthのAPIでも認証部分を手で作る手順は同じなので、応用効きます。
是非色んなAPIを叩いてみて色んなアプリケーションを作ってみてください。
ただし、公開するとAPI Limitをすぐ超えるので、その場合はきちんと認証部分をUIで出来るように作り込みましょう。

Python の Click、結構便利で簡単に CLI ツールが作れて LGTM ですね（← 使い方あってる？）。
その Click を使ってサブコマンドを実装するときに、ファイルを分割して実装するのに少々苦戦したので、そのことについて適当な例を使って実装の流れを書いてみようと思う。
間違いの指摘や、ここもっとこうした方がいいよ等の助言がありましたら、コメント頂けると幸いです。

ちなみに、サブコマンドを複数ファイルに分割して実装しようと思ったのは、「関心事は分離してそれぞれのファイルに分けた方が良いかなぁ」と思ったからである。

TL;DR

subcmd.py

@click.command()
def cmd():
    pass

cmd.py

from subcmd import cmd as subcmd

@click.group()
def cmd():
    pass

cmd.add_command(subcmd)

バージョン

• Python: 3.8.0
• Click: 7.0

サンプルとして実装するもの

良い例が思いつかなかったので、足し算と引き算ができる calculate コマンドを作ってみようと思う。

実行例

$ calculate addition 4 3
7
$ calculate subtraction 4 3
1

ディレクトリ構成

多少順番を変えたり省略したりしているが、下記のような感じである。
ご覧の様(？)に Pipenv を使っている（そろそろ Poetry か何かに乗り換えようかと思っているが、それはまた別のお話）。

calculate
├─ calculate
│  ├─ __init__.py
│  ├─ __main__.py
│  └─ cli
│     ├─ __init__.py
│     ├─ command.py
│     ├─ addition
│     │  ├─ __init__.py
│     │  └─ command.py
│     └─ subtraction
│        ├─ __init__.py
│        └─ command.py
├─ Pipfile
├─ Pipfile.lock
└─ setup.py

下準備

setup.py を書く

Pipenv を好んで使っている理由に、 pipenv install -de . を実行して仮想環境に開発中のツールを編集可能な状態でインストールすることができる点がある（他のパッケージ管理ツールでもできるのかな？）。
そうすることで、仮想環境をアクティベートしたときにコマンドを実行できる。

コマンドを実行するためのエントリーポイントを setup.py に記述する。
下記の例では必要最低限しか記述していない。

setup.py

from setuptools import find_packages, setup

setup(
    name="calculate",
    version="0.0.1",  # 適当
    entry_points={
        "console_scripts": ["calculate=calculate.__main__:main"],
    }
)

パッケージをインストールする

$ pipenv install click
$ pipenv install -de .

コマンドの実装

では、順番に実装を進めていこうと思う。
まずはエントリーポイントを用意する。

__main__.py

def main():
    # TODO: コマンドを呼び出す
    pass

# `python -m calculate` という形で呼び出されたときのことを考慮
if __name__ == "__main__":
    main()

これで、calculate コマンドを実行したときに、__main__.py にある main 関数が呼び出される。
今はまだコマンドを実装していないので、main 関数の中身は空っぽにしておく。

次に main 関数で呼び出す最初のコマンドを実装する。

cli.commnad.py

import click

@click.group()
def cmd():
    pass

@click.group() でデコレートすることで、cmd に対してサブコマンドを追加できるようになる。
calculate コマンドに相当する関数が実装できたので、先ほどの __main__.py の main 関数から呼び出すようにする。

__main__.py

from .cli.command import cmd

def main():
    cmd()

# `python -m calculate` という形で呼び出されたときのことを考慮
if __name__ == "__main__":
    main()

これで実行してみると、さも CLI ツールのような出力が得られる。

$ calculate
Usage: calculate [OPTIONS] COMMAND [ARGS]...

Options:
  --help  Show this message and exit.

やったね?
おっと、まだ目的は達成していなかった。

サブコマンドの実装

CLI ツールとしての基盤(？)が出来上がったので、次に足し算サブコマンドを実装していく。

cli.addition.command.py

import click

@click.command(name="addition")
@click.argument("augend", type=click.INT)
@click.argument("addend", type=click.INT)
def cmd(augend, addend):
    click.echo(augend + addend)

@click.command に name を渡しているが、こうすることでコマンドに任意の名前を付けることができる。
@click.argument では type を渡すことで、バリデーションチェックを行ってくれるようになる。
click.echo で結果を出力している。
ちなみに、augend は足される数、 addend は足す数、という意味である（TDD 本で知った）。

これで、足し算サブコマンドが実装できたので、calculate コマンドに追加してみる。

cli.command.py

import click
from .addition.command import cmd as addition_cmd

@click.group()
def cmd():
    pass

cmd.add_command(addition_cmd)

グループにサブコマンドを追加するには、add_command メソッドを呼び出し、その引数にサブコマンドとして追加する関数を与えるだけである。
ここまで長々と書いてきたが、これだけである。
このメソッドを見つけるまでにとても時間がかかった。。。

引き算サブコマンドの実装も同様である。
コマンドを実装して、

cli.subtraction.command.py

import click

@click.command(name="subtraction")
@click.argument("minuend", type=click.INT)
@click.argument("subtrahend", type=click.INT)
def cmd(minuend, subtrahend):
    click.echo(minuend - subtrahend)

add_command で追加する。

cli.command.py

import click
from .addition.command import cmd as addition_cmd
from .subtraction.command import cmd as subtraction_cmd

@click.group()
def cmd():
    pass

cmd.add_command(addition_cmd)
cmd.add_command(subtraction_cmd)

これで実行してみると、足し算、引き算の結果が得られる。

$ calculate
Usage: calculate [OPTIONS] COMMAND [ARGS]...

Options:
  --help  Show this message and exit.

Commands:
  addition
  subtraction

$ calculate addition 4 3
7

$ calculate subtraction 4 3
1

今度こそ、やったね??

余談

@click.command(name=<another_name>) の name には、記号も与えることができるみたい。
なので、四則演算を実装した場合、次のように書ける。

$ calculate + 4 3
$ calculate - 4 3
$ calculate "*" 4 3
$ calculate / 4 3

* はワイルドカードと認識されるので、クォーテーションで囲まないといけないし、/ は問題なく動くがシェルによってはルートと認識して色がつく。
需要があるかわからないけど、ちょっと面白いなと思ったので載せておきます。
試してないけど、他の記号も色々使えるはず。

まとめ

順番に実装の流れを書いていったので長くなったが、求めていた結果を得るために必要なものは、add_command だった。

間違いの指摘や、ここもっとこうした方がいいよ等の助言がありましたら、コメント頂けると幸いです。

もう数え上げも怖くない ~競プロ数え上げ問題35選~ - Qiitaの記事で紹介されていた問題のうち１つが面白そうだったので解いてみました。ちなみに問題を解きながら考えたことを含めて書いたため、最短で解にたどり着いていません。最短&最良の解を見たい場合は問題文のページから模範解答がリンクされていますのでそちらをご覧ください。

問題

$2N$ 個のマスが左右一列に並んでおり、各マスの色を表す長さ $2N$ の文字列 $S$ が与えられます。
左から $i$ 番目のマスの色は、 $S$ の $i$ 文字目が'B'のとき黒色で'W'のとき白色です。
あなたは異なる $2$ マスを選んで、それらのマスおよびそれらの間にあるマスの色を反転する操作をちょうど $N$ 回行います。 ここで、マスの色を反転するとは、そのマスの色が黒色なら白色に、白色なら黒色にすることです。
ただし、操作を通して同じマスを $2$ 回以上選ぶことはできません。 つまり、各マスがちょうど $1$ 回ずつ選ばれることになります。
$N$ 回の操作終了後に全てのマスを白色にする方法が何通りあるかを $10^9+7$ で割った余りを求めてください。
ここで、条件を満たす $2$ つの方法が異なるとは、$1$ つ目の方法で $i$ 番目に選んだ $2$ つのマスの組と
$2$ つ目の方法で $i$ 番目に選んだ $2$ つのマスの組が異なるような $i (1\leq i \leq N)$ が存在することをいいます。

JSC2019-C Cell Inversion

ちょっとさわってみよう

Bを●または黒、Wを○または白で表します。'反転'には'ひっくり返す'という言葉を使っています。オセロをイメージしてください。
問題文に、●○●●○○○● の例が載っているので、これを解いてみましょう。

手始めに、手順にそってひっくり返してみましょう。

これは ○○○○○○○○ になりませんでした。

もう1回試してみましょう。

これは ○○○○○○○○ になりました。

ひっくり返す手順を、ひっくり返す順に位置をペアにして記載しましょう。
上記の例では、$(2,8)(1,4)(6,7)(3,5)$ になります。このような、ひっくり返す順のペアを並べたものを反転列と呼びましょう。

ひっくり返す順番は関係ない

いくつか試すと、ひっくり返すペアが同じであれば、ひっくり返す順番は最終的に出来上がるマスの模様に関係ないことがわかります。

例えば、$(2,8)(1,4)(6,7)(3,5)$ と $(1,4)(2,8)(3,5)(6,7)$ は同じ結果になります。

なお、以降マスの模様を白黒パターンと呼びます。

反転列を正規化する

反転列はペア同士を入れ替えても最終的に出来上がる模様(白黒パターン)は同じなので、規則を決めて一意に定まるようにします。この操作を反転列の正規化と呼びましょう。

• 小さい方の数字を左に書く
• 左側の数字が小さい順にペアを並べる

$(2,8)(1,4)(6,7)(3,5)$ を上記の規則に則って並べると、$(1,4)(2,8)(3,5)(6,7)$ です。

また、正規化された反転列の左側だけを抜き出したもの(今回だと $\{1,2,3,6\}$ )を左反転列と呼びましょう。

逆に考えてみる

さて、まだどのようにひっくり返せば ●○●●○○○● を全て白にできるのか見当がつきません。

そこで、 全て白のマス( ○○○○○○○○ )を手順にそってひっくり返すとどのようになるのかを、$N$ が小さい範囲で確かめます。

●○●●○○○● を ○○○○○○○○ にすることと、 ○○○○○○○○ を ●○●●○○○● にするのは互いに逆の関係にあります。
○○○○○○○○ を ●○●●○○○● にした後、同じ反転列でひっくり返すと ○○○○○○○○ に戻るからです。

N=1

$N=1$のとき、1通りのひっくり返し方があります。

N=2

$N=2$のときは、3通りのひっくり返し方がありますが、白黒パターンは2種類です。

N=3

$N=3$のときは、15通りのひっくり返し方があり、白黒パターンは5種類です。

規則の発見

ここまででいくつか発見があるので、まとめてみましょう。

最終的に白か黒かは下に引かれている線の数によって決まる

マスの下にある(もしくは自分自身からのびている)線の回数ひっくり返されるので、線の数によって最終的にマスが白になるか黒になるか決まります。
線の数が偶数のとき白、奇数のとき黒になります。

白黒パターンの両端は黒

両端は自身を選択した際にひっくり返され、他のマスに挟まれないので黒です。

クラスタに分かれたとき、各クラスタの両端は黒

クラスタに分かれる、というのは、以下のように全くある区間内で反転が完了している状態です。

クラスタに分かれたとき、各クラスタの両端は黒です。

各パターンの反転列を見ると、左反転列の数字が共通

$N=3$ で、白黒パターンごとの入れ替え列をまとめます。

白黒パターンごとに反転列をまとめると、白黒パターンごとに反転列の左右が同じことに気づきます。例えば、○○○○○○ を ●○○○○● にする反転列には、ペアの左側に$\{1,2,4\}$が現れ、右側に$\{3,5,6\}$が現れます。

重なり合う(L1, R1)(L2, R2)を(L1, R2)(L2, R1)に入れ替えても結果は同じ

下の2つの反転列 $(L_1, R_1)(L_2, R_2)$ と $(L_1, R_2)(L_2, R_1)$ の白黒パターンは同じです。マスの下の線の数が変わらないためです。

互いに線が交差しない反転列が各白黒パターンに1つだけ存在する

線が交差する、とは以下の状況を指します。

ひっくり返すペア同士を線で結んだ時、白黒パターンごとに1つだけ線が交差しない反転列があります。
$N=3$ の場合はそれぞれこの反転列です。

互いに線が交差しない反転列を無交差反転列と呼びましょう。

無交差反転列において、反転のペアは同じ色になる

無交差反転列において、線が繋がっているもの同士は操作終了後同じ色になります。互いに線が交差しない入れ替え列において、ひっくり返すペアは、それ以外のペアに囲まれているか、離れているかのどちらかなので、同時にひっくり返るか同時にひっくり返らないかのどちらかだからです。

隣同士のペアを見つける

無交差反転列では、値が隣同士のペアが少なくとも1つあります。そのペアは白黒パターンでは同じ色になっています。
マスを左から見ていくときに、色が違う場合は隣同士のペアではありません。

スタックに入れてペアを見つけよう

無交差反転列一覧をぼんやり眺めていると、数式のカッコの対応のように見えてきました。

このような対応関係を作成する場合はスタックを使います。

まずは隣同士のペアを見つけましょう。隣同士なので値が連続します。連続した色が見つかるまではスタックに積み上げていきましょう。

●○●●○○○● を例に考えます。

初期状態です。

1番目の●をスタックに入れます。

2番目の○は1番目と違う色なのでスタックに入れます。

3番目の●は2番目と違う色なのでスタックに入れます。

4番目の●はスタックの一番上と同じ色なので、スタックの一番上を取り出し、4番目とペアにします。

5番目の○はスタックの一番上と同じ色なので、スタックの一番上を取り出し、5番目とペアにします。

6番目の○はスタックの一番上と違う色なので、スタックに入れます。

7番目の○はスタックの一番上と同じ色なので、スタックの一番上を取り出し、7番目とペアにします。

最後に、8番目の●はスタックの一番上と同じ色なので、スタックの一番上を取り出し、8番目とペアにします。

ここまでで ●○●●○○○● を ○○○○○○○○ にするには$(1,8)(2,5)(3,4)(6,7)$が1つの解であることがわかりました。実際に確かめて見ると、確かに全て白になります。

例外ケースの対応

上記は反転列が存在するケースですが、ここで例外もみておきます。

スタックの一番下が白

両端は黒でないと、全て白にすることができません。これを言い換えると、スタックに入れる時、一番下は黒でなければならない、となります。一番下に白が入るような入力が与えられた場合、全てのマスを白色にする方法は $0$ 通りです。

最後までマスを読んだときにスタックが空でない

最後までマスを読んだときにスタックが空でない場合は、与えられたパターンになる組み合わせは存在しません。この場合は全てのマスを白色にする方法は $0$ 通りです。

(L1, R1)(L2, R2)を(L1, R2)(L2, R1)に入れ替える場合の数

$(L_1, R_1)(L_2, R_2)$ を $(L_1, R_2)(L_2, R_1)$ に入れ替えても結果は同じです。$(1,8)(2,5)(3,4)(6,7)$ に対して、このような入れ替えが何パターンあるか考えます。

$(1, R_1)(2, R_2)(3,R_3 )(6, R_4)$

の $R_1$ から $R_4$ の部分に$\{4, 5, 7, 8\}$ を入れますが、 $L_i < R_i$ の制約があるので、$R_4$ から考えていきましょう。

$R_4$ : $7$ か $8$
$R_3$ : $4$, $5$, $7$, $8$ のうち1つ (ただし $R_4$ で使用した数字は含まない)
$R_2$ : $4$, $5$, $7$, $8$ (ただし $R_4$, $R_3$ で使用した数字は含まない)
$R_1$ : 残った数字1つ

これを計算すると、$2 \times 3 \times 2 = 12 通り$ です。

ひっくり返す順番の数を掛ける

ここまで、反転列は正規化されていますが、ひっくり返す順番は区別するので、$N!$を掛けます。

$ 12 \times 4! = 288 $

10000000007で割る

最後に、求めた数を $10^9+7$ で割ります。実際のプログラミングでは上記の組み合わせを計算する過程で$10^9+7$で適宜割り、桁があふれないように工夫します。

参考 : 「1000000007 で割ったあまり」の求め方を総特集！ 〜 逆元から離散対数まで 〜 - Qiita

$ 288 \bmod (10^9+7) = 288$

答えの 288 が導かれました。

スタックは無くても解ける

今回はスタックを使用しましたが、スタックが無くても解けます。詳しくは本家の解説をご覧ください。

ソースコード

上記アルゴリズム(を改善したもの)をPython(バージョン3.7)で実装しました。($10^9+7$で割る前まで実装してあります。なので$N$が小さい範囲しか扱えません)

cell_inversion.py

さらなる課題

$N$ 毎に、現れる白黒パターンは、
$N=1$ で $1$種類
$N=2$ で $2$種類
$N=3$ で $5$種類
あることがわかりました。$N=6, 7, 8$と数を大きくしていくとどうなるでしょうか？

あとがき

ざっと書いてみましたが、問題を本格的に考え始めてからアルゴリズムを思いつくまで3日かかりました。(最後ちょっとわからない部分があったので解説を見ました。)本家の解説よりだいぶ遠回りして解いていますが、そのぶん色んな発見があって楽しかったです。
コンテスト参加者はこんな問題をぱっと見で解いてるのかと思うと、末恐ろしいという浅い感想しか出てきません。

背景

現在、競馬×機械学習のために、
ウェブページから競馬情報をスクレイピングすることに取り組んでいます。

入力内容に出走馬の過去の成績を組み込もうとしています。
そこで競馬のレース結果表から出走馬のページのURLを取得、アクセスし、
出走日の過去の成績を取得することにした。

構想

サンプルコード

アイデアは以下の3点。
・取得方法をCSSセレクターで統一
・レース結果表から各馬のURLを取得し、そのURLに掲載されている過去成績を取得
・今回の出走馬は12頭だが、他の馬数でも取得できるように可変長変数を使用

この記事の前に、Python スクレピング 競馬サイトからレース環境を抽出を紹介させていただきました。
こちらはレース環境をBeautifulSoupを使いました。
しかしcssセレクターで統一して処理したほうが気持ちがいいので書き換えました。

scraping_of_race_and_past_horse_result.py

import requests
import lxml.html 
import csv

rlt = [] #結果
horse_URL = []#馬の名前を取得する

#スキャルピングでテキスト取得
def get_scarping_data(key_page,css_select_str,*URL):
    #取得したURLの個数分の情報を取得する
    for i in range(len(URL)):
        #URLから文字列を取得       
        r = requests.get(URL[i])#URLを指定 
        r.encoding = r.apparent_encoding #文字化けを防止
        html = lxml.html.fromstring(r.text) #取得した文字列データ
        # 項目を取得
        for css_id in html.cssselect(css_select_str):        
            #レース結果サイトのレース環境
            if key_page == "condition_in_race_page" :
                #Element番号のテキスト 
                css_id = css_id.text_content()
                #内包表記 抽出する際に"天気"は不変であるので条件指定
                css_id_ = [css_id for t in css_id if "天気" in css_id]
                css_id_ = css_id_[0].split('\xa0/\xa0')
                #リストに行のデータ(リストを追加)
                rlt.append(css_id_)
            #レース結果サイトのレース結果
            if key_page == "result_in_race_page" :
                #Element番号のテキスト 
                css_id = css_id.text_content()
                #改行("\n")をもとに分割する
                css_id = css_id.split("\n")
                #内包表記 空の要素を駆逐する　空文字を除く
                css_id_ = [tag for tag in css_id if tag != '']
                #1位はタイム記録なし #強引に加える必要がある 8番目に0
                if len(css_id_) != 13 : css_id_.insert(8,0)
                #リストに行のデータ(リストを追加)
                rlt.append(css_id_) 
            #レースした馬の過去成績サイト
            if key_page == "horse_race_data" : 
                #Element番号のテキスト 
                css_id = css_id.text_content()
                #取得したElement番号の
                css_id = css_id.split("\n")
                #内包表記 空と"\xa0"と"映像"を除去
                css_id_ = [tag for tag in css_id 
                if tag != '' and tag != "\xa0" and tag != "映像" and tag != "厩舎ｺﾒﾝﾄ" and tag != "備考" ]
                #リストに行のデータ(リストを追加)
                rlt.append(css_id_)        
    #抽出結果
    return rlt

#レースした馬のURLを取得
def get_scarping_past_horse_date(URL):            
    response = requests.get(URL)
    root = lxml.html.fromstring(response.content)
    #1~12までの馬の情報を取得
    for i in range(2,14):#2~13  13 - 2 + 1 = 12 
        css_select_str = "div#race_main tr:nth-child({}) > td:nth-child(4) > a".format(i)
        #レースした馬の情報を取得
        for a in root.cssselect(css_select_str):
            horse_URL.append(a.get('href'))
    #抽出結果
    return horse_URL

#まず本家のレースサイトからレース環境を取得
URL = "https://nar.netkeiba.com/?pid=race&id=p201942100701" 
#レース環境を取得
rlt = get_scarping_data("condition_in_race_page","div#main span",URL)
#レース結果を取得
rlt = get_scarping_data("result_in_race_page","#race_main > div > table > tr",URL)
#レースした馬のURLを取得
horse_URL = get_scarping_past_horse_date(URL)
print(len(horse_URL))#出走馬12に対し、出力結果は12
#取得した馬のURLから過去成績を取得
#項目(1行分)を取得
rlt = get_scarping_data("horse_race_data", "#contents > div.db_main_race.fc > div > table > thead > tr",horse_URL[0])
#項目以外の成績データ
rlt = get_scarping_data("horse_race_data", "#contents > div.db_main_race.fc > div > table > tbody > tr",*horse_URL)
#CSVファイルに保存
with open("scraping_of_race_and_past_horse_result.csv", 'w', newline='') as f: 
    wrt = csv.writer(f) 
    wrt.writerows(rlt) #抽出結果の書き込み

実行結果

反省点

・レース結果から馬の情報が掲載されているURLを取得し、過去成績も取得できた。
・過去成績を取得できたが最新の情報も掲載されているので、レースした日を基準に遡った成績を取得する必要がある。
・取得できた情報をニューラルネットワークに入力できるレベルまで行っていない。
・レースした日を自動で取得する。毎回調べない
→方法として、レース情報数値(URLの末端の12桁)を規則性を考慮してインクリメント。

CSSセレクタの調べ方

スクレイピングするうえで、以下のコード(文字列)を使いこなす必要がある。
(#contents > div.db_main_race.fc > div > table > tbody > tr')
↑ここの部分をどのようにして取得するのかをまとめる。

how_to_search_css_sector.py

for h in html.cssselect('#contents > div.db_main_race.fc > div > table > tbody > tr'):#スクレイピング箇所をCSSセレクタで指定

初めはただコピペしたり、そっれっぽいコードを書いてトライしたが、うまく行かない。。。
調べてみるとCSSが絡んでいる。
欲しい項目はどのCSSセレクターで構成されたいるかを知るために、
ChromeのCopy Css Selectorのツールを使う。

インストールすると右クリックでCopy Css Selectorの項目ができる。
欲しい項目でCopy Css Selectorを実行し、テキストに張り付けて確認するとCSSセレクターを知ることができる。
※うまくいかないこともあるかもしれない。

他にもディベロッパー ツールを用いて、
欲しい情報範囲で　Copy → Copy Selector　でCCSセクターをコピーできる。
たまにうまくいきます。

この2つを駆使すれば、SCCセレクターを取得できるはずです。

記事の内容

子供の頃に台風ってこんなに頻繁に発生してたっけ？とか10月ぐらいに台風って日本に来てたっけ？という漠然とした疑問がありました。
気象庁が公開している台風に関するデータを使ってこのあたりの疑問を解決出来たらないいなと思い、色々やってみました。

データ取得元

政府が運営しているデータカタログサイトから取得しました。

データカタログサイト

気象予報_天気予報・台風の資料
このページの「台風位置表のCSVデータ」を使用します。

この記事ではリンク先に掲載されている2001年から2019年までのCSVデータを使用しました。

foliumを使って台風の経路を地図上に書いてみる

準備

CSVデータのヘッダーが日本語だったので英語に変えました。だいぶ適当です

year,month,day,hour(UTC),typhoon_no,typhoon_name,rank,latitude,longitude,central_pressure,max_wind_speed,50KT_LDD,50KT_LD,50KT_MA,30KT_LDD,30KT_LD,30KT_MA,landing

コード

あまり、工夫したところはありません

typhoon_plot.py

import folium
import pandas as pd

pd.options.display.precision = 3

# 開始と終了の年
start_year = 2001
end_year = 2019

# 1年単位で処理する
for year in range(start_year, end_year+1):
    print('# ' + str(year) + ' year start')
    # データの読み込み
    df = pd.read_csv('./typhoon/table' + str(year) + '.csv',encoding="SHIFT-JIS")
    typhoon_names = list(df['typhoon_name'])
    latitude = list(df["latitude"])
    longitude = list(df["longitude"])
    landing = list(df["landing"])

    map = folium.Map(location=[35.6611098,139.6953576], zoom_start=3)

    target_typhoon = typhoon_names[0]
    target_location = []
    color = 'blue'
    typhoon_count = 0

    # CSVデータを1行ずつ処理する
    for lt, lo, name, land in zip(latitude, longitude, typhoon_names, landing):
        # 台風の名前が変わったら処理中の台風の情報を地図に書く
        if name != target_typhoon:
            map.add_child(folium.PolyLine(locations=target_location, color=color))
            target_location = [[float(lt), float(lo)]]
            target_typhoon = name
            typhoon_count += 1
        else:
            target_location.append([float(lt), float(lo)])

    map.add_child(folium.PolyLine(locations=target_location, color=color))
    map.save('./output/' + str(year) + '_typhoon_location.html')
    print('# Number of typhoon : ' + str(typhoon_count))
    print('# ' + str(year) + ' year end')

出力結果

出力結果は年毎に「output」フォルダの配下に出力しています。
出力したHTMLを開くとこんな感じになります。

出力結果の画像を2001年から2019年まで表示するGIFを作ってみました。

日本に上陸した台風は線の色を変えてみる

やはり、日本に上陸した台風は目立たせたいですよね。
データの最終項目に「上陸」という項目がありました。上陸していたら１、していなかったら０です。
このデータを使って色分けをします。

注意しないといけないのが、この「上陸」は「日本に」ではありません。
そのため、日本かどうかの判定を入れる必要がありました。

コード

typhoon_plot.py

import folium
import pandas as pd

pd.options.display.precision = 3

# 開始と終了の年
start_year = 2001
end_year = 2019

# 日本の東西南北の端の座標
east_end = 153.5911
west_end = 122.5601
north_end = 45.3326
south_end = 20.2531

def is_japan_randing(lt, lo, randing):
    if south_end <= lt and lt <= north_end and west_end <= lo and lo <= east_end and randing == 1:
        return True
    else:
        return False

# 1年単位で処理する
for year in range(start_year, end_year+1):
    print('# ' + str(year) + ' year start')
    # データの読み込み
    df = pd.read_csv('./typhoon/table' + str(year) + '.csv',encoding="SHIFT-JIS")
    typhoon_names = list(df['typhoon_name'])
    latitude = list(df["latitude"])
    longitude = list(df["longitude"])
    landing = list(df["landing"])

    map = folium.Map(location=[35.6611098,139.6953576], zoom_start=3)

    target_typhoon = typhoon_names[0]
    target_location = []
    typhoon_count = 0
    color = 'blue'

    # CSVデータを1行ずつ処理する
    for lt, lo, name, land in zip(latitude, longitude, typhoon_names, landing):
        # 台風の名前が変わったら処理中の台風の情報を地図に書く
        if name != target_typhoon:
            map.add_child(folium.PolyLine(locations=target_location, color=color))
            target_location = [[float(lt), float(lo)]]
            target_typhoon = name
            typhoon_count += 1
            color = 'blue'
        else:
            target_location.append([float(lt), float(lo)])
            if is_japan_randing(lt, lo, land):
                color = 'red'

    map.add_child(folium.PolyLine(locations=target_location, color=color))
    map.save('./output/' + str(year) + '_typhoon_location.html')
    print('# Number of typhoon : ' + str(typhoon_count))
    print('# ' + str(year) + ' year end')

「is_japan_randing」の中で日本国内の判定をしています。
凄く雑な気がしますが、緯度が日本の最北端、最南端の範囲か？経度が最西端、最東端の範囲か？を判定しています。

結果

日本上陸している線だけが赤くなりました。
九州や本州にかなり近づいている線がありますね。以下、お天気.comの引用です。

台風の上陸とは台風の中心が北海道、本州、九州、四国の海岸に達した場合を言います。 （沖縄は台風の上陸とは言わず、通過と言う。） 台風の接近とは半径３００ｋｍ以内に入る事を言います。

こう考えると、この線は中心の座標なんだということが改めて分かりますね。

終わりに

近年こんなに台風って日本に上陸してたっけ？という疑問からこれをやってみましたが、上陸数とかは2001年以降だとあまり変わってないようですね。
次は台風の規模も可視化してみようと思います。

Databricks とは、大量のデータを並列に処理するアプリを作るためのサービスだ。Apache Spark の開発者によって Apache Spark のマネージド・サービスとして開発された。ここ数日勉強しているが特に引っかかった点を書く。

一言でいうと Apache Spark とは、表形式のデータを加工するコードを自動的に並列処理に変換して並列で実行してくれるという物だ。開発者があたかも Jupyter 上の Pandas でコードを書いているような気分で巨大データを並列に加工出来る。機械学習関のライブラリも含まれているので、データの前処理から分析や予測まですべて Spark で行う事が出来る。

また、前身の(?) Hadoop との違いは、ストレージを切り離した所にある。Spark はデータ処理に特化してする事で様々な他のストレージと組み合わせて使えるようになった。

主な情報源

• Databricks ドキュメント: https://docs.databricks.com/
• 参考書: Spark: The Definitive Guide http://shop.oreilly.com/product/0636920034957.do
  • Apache Spark についての本だが、勉強のお勧めは Databricks Community Edition を使う事とある
    • https://community.cloud.databricks.com/

以下に Databricks で分かりにくい事があったので記録する。

CLI

分かりにくくは無いが CLI が無いとかなり面倒なのでインストールは必須。インストールには pip3 を使う。databricks configure で接続情報を設定すると使えるようになる。

# 設定
$ pip3 install databricks-cli
$ databricks configure --token
Databricks Host (should begin with https://): https://hogehoge.cloud.databricks.com/
Token: (GUI で作った token を入力)

# 動作確認
$ databricks fs ls

databricks コマンドには様々なサブコマンドがあるが、databricks fs には最初から短縮形として dbfs が用意されている。

Secret

パスワードなどは scope ごとに Secret に保存する。

scope 作成

databricks secrets create-scope --scope astro_snowflake

scope 表示

databricks secrets list-scopes

scope に secrets 追加

databricks secrets put --scope hoge_scope --key User
databricks secrets put --scope hoge_scope --key Password

secrets 表示

databricks secrets list --scope astro_snowflake

DBFS と Workspace の２つのファイル置き場

Databricks で利用するファイルは以下の二箇所がある。

• Databricks File System (DBFS)
  • Python プログラムからは /dbfs にマウントされているように見える。
  • プログラムで利用するファイルはここに置く。
  • cli からは databricks fs コマンドでアクセス出来る。
  • 特に、/FileStore/ 以下はファイルのアップロードなど UI からアクセスするのに使う。
• Workspace
  • Notebook はここに置く。
  • cli からは databricks workspace コマンドでアクセス出来る。

なんと、プログラムから直接 Workspace にアクセスする事は出来ない！

SQL Table と DataFrame の関係

Databricks では、表にアクセスするのに Python や SQL を使う事が出来る。同じ表を別の言語で使うには表を SQL から見える所に Table として登録する必要がある。

Table にはどこからでもアクセス出来る global と、同じ Notebook からだけアクセス出来る local がある。

Python の DataFrame を "temp_table_name" という local table として登録する。これで SQL から参照出来る。

df.createOrReplaceTempView("temp_table_name")

Python の DataFrame を "global_table_name" という global table として登録する。これで SQL から参照出来る。global table は Web UI の Data から参照出来る。

df.write.format("parquet").saveAsTable("global_table_name")

"temp_table_name" という名前で登録されてある table を Python の DataFrame として読む。

temp_table = spark.table("temp_table_name")

global table は DBFS 上に保存される。保存される場所は DESCRIBE DETAIL で出てくる location で確認出来る。

DESCRIBE DETAIL `global_table_name`

例えば dbfs:/user/hive/warehouse/global_table_name

Notebook のマジックコマンド

https://docs.databricks.com/notebooks/notebooks-use.html

マジックコマンド

• %fs dbutils.fs を簡単に使う
  • 例: %fs ls /databricks-datasets/structured-streaming/events/
  • dbutils とは: https://docs.databricks.com/data/databricks-file-system.html#dbfs-dbutils
• %sh shell コマンド
  • 例: %sh head -n 5 /dbfs/databricks-datasets/structured-streaming/events/file-0.json
• %md: Markdown

はじめに

今話題のオープンデータですが、必ずしもCSVなどの生データではない事で、微妙に扱いにくいデータに仕上がっている事があります。
保有しているデータを新たに出している事自体は当然褒められるべき事で、「ほんとだったら生データ欲しいなァ〜」くらいの気持ちで期待している訳ですね。
そんな訳で、PDFで出したから悪いわけではありません、が、データの活用に際しては、PDFではなくより機械で読みやすい形式にする必要があります。
今回は、札幌市交通局が公開した朝ラッシュ時間帯の車内混雑状況についてのPDFデータを題材に、CSVデータを錬成する手順を紹介してみます。

方針

上記ウェブサイトで得られるPDFファイルはこのような書式です

一応データは構造化されていて、エクセルをPDF出力した感じです。
なので「構造から解析路線」が思い浮かびます、が、表に格納されている日本語が解析時に化けてしまうなど手間がかかりそうでした。
そんな訳で次に、以下の画像のような手順を考案。

画像処理ソフト使うのが面倒だったのでiPadで手書きです。
要は、いくつかPDFを見ると、左表も右表も、左上のセルは常に同じ位置にあって、セルのサイズも全部同じという事がわかった訳です。
ならば、画像で言う赤点の箇所のRGB値（＝色）を抽出・混雑度データに変換してやれば、望みのデータとなる訳ですね。

結論

こんな風に実装できました
https://github.com/Kanahiro/sapporo_metro_analyze/

こんな風にCSVを出力します

赤いとこは4、白は0で青は2になっている事がわかりますね。
では、PDF読み込み→画像に変換→指定ピクセルの色を取得→色から当該セルの混雑度データを生成→CSVファイルに、という手順を追ってみます。

PDFから画像に変換

参考：PythonのPDF処理まとめ(結合・分割, 画像変換, パスワード解除)

pdf2imageを利用します。使い方は上記記事を参照。
なお、記事ではpip install popplerとありますが、現在はpipではインストール出来ません。
Linuxなら以下です。他のOSの説明は省略します。

sudo apt install poppler-utils

画像から指定ピクセルのRGB値を取得

pdf2imageで読み込んだデータは、numpyのarrayに変換する事が出来ます。
つまり、ピクセル構造と一致した2次元配列に、RGB配列を突っ込んで、３次元配列となります。

#convert_from_pathはpdf2imageの関数
pdf_images = convert_from_path(pdffile)
img_array = np.asarray(pdf_images[0])
'''
img_arrayのサンプル

[[[255 255 255]
  [255 255 255]
  [255 255 255]
  ...
  [255 255 255]
  [255 255 255]
  [255 255 255]]

 [[255 255 255]
  [255 255 255]
  [255 255 255]
  ...
  [255 255 255]
  [255 255 255]
  [255 255 255]]

 [[255 255 255]
  [255 255 255]
  [255 255 255]
  ...
  [255 255 255]
  [255 255 255]
  [255 255 255]]

（中略）

 ...

 [[255 255 255]
  [255 255 255]
  [255 255 255]
  ...
  [255 255 255]
  [255 255 255]
  [255 255 255]]

 [[255 255 255]
  [255 255 255]
  [255 255 255]
  ...
  [255 255 255]
  [255 255 255]
  [255 255 255]]

 [[255 255 255]
  [255 255 255]
  [255 255 255]
  ...
  [255 255 255]
  [255 255 255]
  [255 255 255]]]

PDFの端っこは白いので当たり前ですが[255 255 255]すなわち白ばかりですね。
ピクセル単位に配列に格納出来ました。

特定ピクセル(data_of_pixel)へのアクセスは

x = START_CELL[0] + c * CELL_SIZE[0] #x座標
y = START_CELL[1] + r * CELL_SIZE[1] #y座標
data_of_pixel = img_array[y][x]

となります。これで、PDFを変換した画像の、特定の箇所のRGB値を取得出来ました。

混雑度の判定

上記の画像から、白、水色、青、黄色、赤の順で混雑度が高い事を示す事がわかります。
ところが、右上の凡例とデータエリアで、若干RGB値が異なるPDFに仕上がっていました。
幸い、段階ごとにそれなりに色の違いが大きいので、ここは凡例とデータエリアのセルとのRGB値の差の大きさで判定したいと思います。

#混雑度凡例のRGB値
CROWD_RGBs = [
    [255, 255, 255],
    [112, 200, 241],
    [57, 83, 164],
    [246, 235, 20],
    [237, 32, 36]
]

def rgb_to_type(rgb_list)->int:
    #色差の閾値
    threshold = 50
    color_array = np.asarray(rgb_list)
    for i in range(len(CROWD_RGBs)):
        crowd_rgb_array = np.asarray(CROWD_RGBs[i])
        color_dist = abs(color_array - crowd_rgb_array)
        sum_dist = color_dist.sum()
        if sum_dist < threshold:
            return i #0 - 4 混み具合

この関数にRGB値のlistを渡すと、混雑度を0-4の整数値で返します。
何をしているかというと、先ほど取得したピクセルのRGB値と、混雑度凡例のRGB値を比較しています。
RGB値それぞれの差の絶対値の総和を色差と定義し、その差が50以内の場合、混雑度を確定しています。

これで、すべてのセルの混雑度を判定し、CSVにすれば冒頭のCSVデータが完成します。

終わりに

こういった処理はオープンデータ界隈ならあるあるネタなんでしょう。「餅から米をつくる」とおっしゃっていた方もいましたが、これで私も錬米術師を名乗ってもよいでしょうか？ただやはり、一次データ（餅）がない限り米の錬成は出来ない訳なんで、そこには感謝しかない訳です、いつもありがとうございます。データの量は増えてきたので、次はそのデータの質の段階になっていくのかな…。

はじめに

メタプログラミングとは

メタプログラミング (metaprogramming) とはプログラミング技法の一種で、ロジックを直接コーディングするのではなく、あるパターンをもったロジックを生成する高位ロジックによってプログラミングを行う方法

出展 : wikipedia

要するに、コードを生成するコードを書いて複雑なコードを書いちゃいましょう！ってこと（だと私は解釈しています。間違っていたら教えてください）。

概要

Pythonの記号計算用ライブラリであるSympyを使って数式を生成し、C言語にコンバートします。SympyはMathematicaのように数式処理に使われることが多いようですが、数式を様々な言語にコンバートする機能も兼ね備えています。私はあまり使ったことがありませんが、ヘッダファイル付きで関数を出力する機能もあります。

この記事では入門的に、n×n行列の掛け算のコードをSympyでメタプログラミングしていきます。
Sympyの記法や使い方は以下にまとめられているので、興味のある方は読んでみてください。

• Welcome to SymPy’s documentation!
• 3.2. Sympy : Python での代数計算 （日本語）

インストール

Sympyのインストールはpipでできます。

pip install sympy

ダウンロードの場合は以下から。

• https://www.sympy.org/en/download.html

インポート

まずはSympyをインポートします。

from sympy import *

今回は行列計算を実装したいので、Matrixもインポートしておきます。

from sympy import Matrix

これでMatrixも宣言できるようになりました。

シンボルの作成

数式に用いる変数はsymbolとして事前に宣言しておく必要があります。例えば、xとyをsymbolとして宣言する場合は

x = symbols('x')
y = symbols('y')

とすればOKです。今回は行列計算のために、n×nのMatrixの要素としてsymbolを定義していきます。

n = 3
A = Matrix([[symbols('a['+str(j)+']['+str(i)+']') for i in range(n)] for j in range(n)])
B = Matrix([[symbols('b['+str(j)+']['+str(i)+']') for i in range(n)] for j in range(n)])

これで行列A,Bがそれぞれ以下のように宣言されました。ここではn=3としました。

⎡a[0][0]  a[0][1]  a[0][2]⎤
⎢                         ⎥
⎢a[1][0]  a[1][1]  a[1][2]⎥
⎢                         ⎥
⎣a[2][0]  a[2][1]  a[2][2]⎦

⎡b[0][0]  b[0][1]  b[0][2]⎤
⎢                         ⎥
⎢b[1][0]  b[1][1]  b[1][2]⎥
⎢                         ⎥
⎣b[2][0]  b[2][1]  b[2][2]⎦

行列の出力はpprint(A)と記述するとこのように視覚的にわかりやすい形になります。

行列の掛け算

先ほどの行列A,Bの積を行列Cとします。行列の掛け算は簡単に、以下でOKです。

C = A*B

今度はprint(C)で中身を見てみると

Matrix([[a[0][0]*b[0][0] + a[0][1]*b[1][0] + a[0][2]*b[2][0], a[0][0]*b[0][1] + a[0][1]*b[1][1] + a[0][2]*b[2][1], a[0][0]*b[0][2] + a[0][1]*b[1][2] + a[0][2]*b[2][2]], [a[1][0]*b[0][0] + a[1][1]*b[1][0] + a[1][2]*b[2][0], a[1][0]*b[0][1] + a[1][1]*b[1][1] + a[1][2]*b[2][1], a[1][0]*b[0][2] + a[1][1]*b[1][2] + a[1][2]*b[2][2]], [a[2][0]*b[0][0] + a[2][1]*b[1][0] + a[2][2]*b[2][0], a[2][0]*b[0][1] + a[2][1]*b[1][1] + a[2][2]*b[2][1], a[2][0]*b[0][2] + a[2][1]*b[1][2] + a[2][2]*b[2][2]]])

と、きちんと計算出来ていることがわかります。

C言語へのコンバート

C言語へのコンバートは

ccode(<数式>, <代入する変数>, standard='C99')

と記述します。今回は行列Cの全要素をコンバート・出力したいので、

for i in range(n):
    for j in range(n):
        idx = i*n+j
        code = ccode(C[idx],assign_to=('c['+str(i)+']['+str(j)+']'), standard='C89')
        print(code)

としました。idx = i*n+jは、i, jを一次元に変換する式です。上記コードを実行すると以下が出力されます。

c[0][0] = a[0][0]*b[0][0] + a[0][1]*b[1][0] + a[0][2]*b[2][0];
c[0][1] = a[0][0]*b[0][1] + a[0][1]*b[1][1] + a[0][2]*b[2][1];
c[0][2] = a[0][0]*b[0][2] + a[0][1]*b[1][2] + a[0][2]*b[2][2];
c[1][0] = a[1][0]*b[0][0] + a[1][1]*b[1][0] + a[1][2]*b[2][0];
c[1][1] = a[1][0]*b[0][1] + a[1][1]*b[1][1] + a[1][2]*b[2][1];
c[1][2] = a[1][0]*b[0][2] + a[1][1]*b[1][2] + a[1][2]*b[2][2];
c[2][0] = a[2][0]*b[0][0] + a[2][1]*b[1][0] + a[2][2]*b[2][0];
c[2][1] = a[2][0]*b[0][1] + a[2][1]*b[1][1] + a[2][2]*b[2][1];
c[2][2] = a[2][0]*b[0][2] + a[2][1]*b[1][2] + a[2][2]*b[2][2];

行列Cの各要素が計算され、C言語の形で出力されました。ちゃんと ; も付いていますね。あとはこれを自分のコードに張り付ければOKです。

ソースコード

今回作成した行列計算コードです。

from sympy import *
from sympy import Matrix
from sympy.utilities.codegen import codegen
n = 3
A = Matrix([[symbols('a['+str(j)+']['+str(i)+']') for i in range(n)] for j in range(n)])
B = Matrix([[symbols('b['+str(j)+']['+str(i)+']') for i in range(n)] for j in range(n)])
pprint(A)
pprint(B)
C = A*B
print(C)
for i in range(n):
    for j in range(n):
        idx = i*n+j
        code = ccode(C[idx],assign_to=('c['+str(i)+']['+str(j)+']'), standard='C89')
        print(code)

まとめ

Sympyを使って行列の掛け算を計算するコードをメタプログラミングしました。このくらいの計算ならメタプログラミングするまでもなく書けますが、「元の式をtで何階偏微分して～」みたいな複雑な数式を実装する時には非常に有用です。（行列掛け算もfor文がなくなった分少し実行速度が上がるメリットはある？）

ご参考になれば幸いです!

どういう記事か？

• 最近「ゼロから作るニューラルネットワーク」という本を読んで学んだニューラルネットワークの4つの最適化手法についてまとめようと思います
• 基本的には本がベースになっていますが、どうイメージしたらわかりやすいか？ということに重点を置いて個人的な解釈をまとめた記事です。なので、厳密な数式や応用などは載せていません。
• ネットで調べるとこの記事よりももっと詳しい説明が腐るほど出てきますが、一方で本とほぼ同じ説明をそのまんまQiitaにまとめた記事も散見されたので、自分なりの解釈をまとめておく意味も少しはあるかなと思い、記事を書くに至りました
• Momentumという手法の説明にかなり偏ってしまいましたが、ご容赦ください。。
• 僕はニューラルネットワークについては初学者、初心者なので、間違っている箇所やアドバイスなどがあればコメントもらえると喜びます

ゼロから作るニューラルネットワーク

今や機械学習の専門でない一般の人でも名前だけは知っている「ニューラルネットワーク」「ディープラーニング」。
機械学習の１つで、画像認識などによく使われているやつですね（Google画像検索みたいなやつ）
このニューラルネットワークを勉強する際におそらく最も日本で読まれている本が「ゼロから作るDeep Learning」です。
よく読まれているだけあって、入門書としては非常に良書だと思います。（何様なんだ）

機械学習における"学習"

ニューラルネットワークに限らず、機械学習の目的は学習データ（教師データ）をもとに学習し、未知のデータの予測を行う、というものです。
例えば画像認識で言えば猫の画像を大量に学習させて、新しく画像を見せた時に「その画像か猫かどうか？」を答えさせる、というものですね。
通常の機械学習であれば「特徴量を指定して学習しなければならない」 という制約がありますが、ディープラーニングの場合は学習する中で画像から勝手に特徴量を見つけてくれるという利点があります。

最適化ってなーに？

学習とはいったい何をしているのか？というと、最適な重みパラメータを探す作業に当たります。
最適なパラメータが求まればそれを元に未知のデータから正解を予測できるようになります。
これを最適化(optimize)と言います。
最適化のためには損失関数（コスト関数）と呼ばれる関数を定義し、最小になるパラメータを探します。
損失関数としては最小二乗誤差やLogLoss (交差エントロピー誤差)などが用いられ、いずれも小さい値を取るほどうまく分類（or 回帰）できているということになります。
ではどうやって損失関数の最小の場所を探すか？を考えます。
この本では以下の4種類の最適化手法が紹介されています。

• SGD（勾配降下法）
• Momentum
• AdaGrad
• Adam

SGD

勾配降下法（Gradient Stochastic Descent）の略。
最適化手法の中で一番スタンダードなやつで、以下の式で重みを更新します

\boldsymbol{W} \leftarrow \boldsymbol{W} - \eta\frac{\partial{L}}{\partial{\boldsymbol{W}}}

損失関数の傾きを計算して損失関数が小さくなるように更新していきます。
「坂道を下に下っていけばいつかは一番低いところに辿りつけるだろう」という考えです。
$\eta$は学習率と呼ばれるやつで、小さすぎるとなかなかたどり着かず、大きすぎると一度に大きく動きすぎてとんでもない方向に飛んで行ってしまったりします。
なので、適切な値を最初に設定してあげる必要があります。
個人的にはゴルフのようなものだと思っていて、ちょこっとの力でしか打たないといつまでたってもホールにたどり着かないし、かといって力任せに打ちすぎるとホールを通り越してとんでもない方向にいってしまいます。
そんなSGDですが、常にうまくいくとは限りません。例えば$f(x,y)=\frac{1}{20}x^{2}+y^{2}$のような関数を考えます。お椀型をx軸方向にビヨっと引き伸ばしたような形ですね。
この関数でSGDをやってみた結果が以下です。

縦方向の勾配が急なので縦にギザギザしながら原点（損失関数が最小の点）に向かっていってるのがわかります。かなり効率の悪い動きをしていますね。これを改善したのが以下に出てくる３つの方法です。

Momentum

Momentumとは物理学で出てくる運動量のことです。まずは重みの更新式を見てみましょう。

 \boldsymbol{v} \leftarrow \alpha\boldsymbol{v} - \eta\frac{\partial{L}}{\partial{\boldsymbol{W}}}\\
\boldsymbol{W} \leftarrow \boldsymbol{W} +\boldsymbol{v}

式を見ただけだとピンと来ないかもしれないので、まずは先ほどの関数に適用した結果を見てみましょう。

滑らかな動きになりましたね。なんとなくお椀の中をボールが転がったような軌跡を描いていて、自然です。
なぜ滑らかな動きになるのか？以下に簡単な概念図を書いてみました。

一つの例を考えてみましょう。①が最初の学習、②が2回目の学習を表すとします。
最初に$\frac{\partial{L}}{\partial{\boldsymbol{W}}}$を計算した結果、①のように上方向のベクトルとなったとします。その方向にしたがって$\eta$の大きさを掛けた分だけ移動したとします。
ここで②の重みを更新を考えます。再度Lの微分を計算した結果、今度は斜面の下り坂の方向が上ではなく右方向に変化していたとしましょう。SGDの場合だと①と同じように今度は②で右方向に進みます。
結果として直角に右に移動していますね。これがSGDのガタガタの原因で、不自然は軌道になってしまっています。
ここでMomentumの式を確認してみましょう。更新式の中に$\alpha{v}$の項がありますね。
これは①の時の$\frac{\partial{L}}{\partial{\boldsymbol{W}}}$にあたります。
簡単のため$\alpha=1$だったとしましょう。
すると$\alpha\boldsymbol{v}-\eta\frac{\partial{L}}{\partial{\boldsymbol{W}}}$の部分は①と②のベクトルの和になっていることがわかります。
すると②の方向は右ではなく斜め45度で右上に進むことになります。
このように、前回の重み更新時の傾きを次回の重み更新に反映させることで、動きが滑らかになります。
αは前回の更新時の傾きをどれくらい影響させるかというパラメータになります。
αが0なら前回の重み更新時の傾きは一切反映しないので右方向、つまりSGDと同じになります。

次に、なぜこの手法がMomentumと呼ばれるのかを考えてみます。物理学（ニュートン力学）では運動量は質量(m)と速度のベクトル($\boldsymbol{v}$)の積になります。
運動量と力積（$F\Delta{t}$）の関係という以下のような関係式が成り立ちます。

 m\boldsymbol{v_{before}}+F\Delta{t} = m\boldsymbol{v_{after}}

これも図で考えてみましょう

上方向に運動量$m\boldsymbol{v_{before}}$で運動している物体が右方向に力積を受けた場合に物体の運動がどう変化するかを表します。
もし力積の大きさが最初の運動量と同じであれば斜め45度に右上に、力積の大きさがそれより大きくなると、より右側に物体は弾き飛ばされます。
これはさっきの概念図と非常に似ているということがわかると思います。
Momentumによって描かれた軌跡のことを先ほどお椀の中をボールが転がるようにと表現しましたが、実際の力学の物理法則と似たような式になっているためこのような軌跡になるということがわかると思います。
イメージ的には上に進んでいるところに急に右方向に力がかかっても急に直角には曲がれないよ！というのを重みの更新に適用したという感じですかね。

AdaGrad

Momentumの説明が非常に長くなってしまいましたが、残りの方法についてもさらっと説明したいと思います。
SGDの説明の時に、重みの更新をゴルフに例えました。SGDでは毎回の学習で同じ学習率で更新します。つまりゴルフで言えば毎回同じ強さで球を打っていることになります。
しかし、ゴルフの場合、一打目は強く打って飛距離を伸ばし、ホールが近づいてきたらだんだん弱く打って微調整していきますよね？
同様のことがニューラルネットワークでも適用可能です。つまり最初の学習時は大きく重みを更新し、だんだんと更新する大きさを減らしていくという方法です。
この手法をAdaGradといいます。
以下の式をみてみましょう。

 \boldsymbol{h} \leftarrow \boldsymbol{h} + \left(\frac{\partial{L}}{\partial{\boldsymbol{W}}}\right)^{2}\\
\boldsymbol{W} \leftarrow \boldsymbol{W} -\frac{1}{\sqrt{\boldsymbol{h}}}\eta\frac{\partial{L}}{\partial{\boldsymbol{W}}}

下の式で$\sqrt{\boldsymbol{h}}$で割っていますね。これは毎回重みを更新するごとに更新する大きさが小さくなっていくことを意味します。
上の式の$\left(\frac{\partial{L}}{\partial{\boldsymbol{W}}}\right)^{2}$は行列の全要素を２乗したものを意味します。（内積ではないので、スカラーにはならない）
つまり、傾きが急であるほど大きく重みを更新し、だんだん小さくしていくということですね。
非常に理にかなっています。
これを実際に先ほどの関数に適用したのが以下になります。

今回の例だと非常にうまくいっているように見えますね。

Adam

最後にAdamについて紹介します。
Adamは簡単に言えばMomentumとAdaGradのいいこと取りしたやつです。
数式含め、ここでは詳しく触れません。
例によって関数に当てはめると以下のような図になります。

これもSGDと比べるといい感じに最小値にたどり着けてますね。

まとめ

• ニューラルネットワークで用いる４つの最適化手法を（かなりざっくり）紹介しました
• 「結局どれ使いばいいの？」という疑問ですが、万能な方法は存在しないので臨機応変に使い分ける必要がありそうです
• SGDは今でもよく用いられていますが、最近ではAdamが人気なようです

備考

• 今回図をplotするのに使ったコードは以下のgitで公開されています
• https://github.com/oreilly-japan/deep-learning-from-scratch

Pythonとは？

• Pythonの特徴

動的な型付け	型宣言をプログラム実行前にしなくて良い
ガベージコレクション	不要なメモリは自動的に解放
マルチパラダイム	オブジェクト指向、命令型、手続き型、関数型など複数の特性を併せ持つプログラミング言語
外部ライブラリが豊富	機械学習用の「TensorFlow」「Chainer」や Web開発用の「Django」「Flask」など様々
読みやすい	簡単な短いコードによって書けて読みやすい
多くのプラットフォームに対応	Windows,Mac,Linux,ios,Android...

C言語でのプログラム実行まで

たった１行の出力でもいろいろしなくちゃいけない...
＊もちろんC言語も高速であったり自由度が高かったりと優れているところがたくさんあります。

Hello.c

//ライブラリの読み込み
#include <stdio.h>
//main関数の定義
int main(int argc, char *args[])
{
    printf("Hello, world!\n");
    return 0;
}

$  gcc -o hello hello.c　　　　　　　//コンパイル
$  ./hello　　　　　　　　　　　　　　 //実行
Hello, world!
$

Pythonでのプログラム実行まで

たったこれだけで動きます。
Pythonはインタプリタ言語といって
コンピュータがプログラミング言語を機械語に翻訳する作業を
プログラム実行時に同時作業で行ってくれます。
これにより、機械語への変換を一括で行うコンパイルという作業の指示をする必要がありません。

Hello.py

print("Hello, world!")

$ python Hello
Hello, world!

---

実際に書いてみよう！

Google Colaboratory
* ファイルをクリック

* ノートブックを新規作成をクリック

* 以下のコードを記入

print("Hello, world!")

• ボタンを押す

変数を使ってみよう

msg = "Hello, world!"
num = 100
ADMIN_EMAIL = "hello@mail.com"
squares = [1,4,9,16,25]
print(msg)
print(num)
print(squares)
print(ADMIN_EMAIL)


msg = "I'm Pretender"
num = num + 50
squares += [36, 49, 64, 81, 100] 
print(msg)
print(num)
print(squares)

Hello, world!
#numの出力
100
#squaresの出力
[1, 4, 9, 16, 25]
#ADMIN_MAILの出力
hello@mail.com
#更新後のmsgの出力
I'm Pretender
#更新後のmsgの出力
150
#更新後のsquaresの出力
[1, 4, 9, 16, 25, 36, 49, 64, 81, 100]

数値の四則演算など

x = int(input("x="))
print(x+3)
print(x-3)
print(x*3)
print(x/3)
msg += "GoodBye!"
print(msg)

文字列に値を埋め込む

name = "Iwana"
height = 173
#やり方１　文字列フォーマット指定
print("%s くんの身長は、 %d cmです。" %(name,height)) 
#やり方２ {}で指定
print("{1} くんの身長は、 {2} cmです。".format(name,height))

フォーマット指定子	型
%d	整数
%f	浮動小数
%s	文字列
%o	8進数
%x	16進数

フォーマット指定子のサイト

課題

名前とBMIを計算して表示するプログラムを書いてみよう。

条件分岐をしてみよう

if文

score=int(input("score ? "))
if score >= 90:
    print("A+")
elif score >= 80:
    print("A")
elif score >= 70:
    print("B")
elif score >= 60:
    print("C")
else:
    print("You're failed.")

Pythonに限らずプログラミング言語では基本的に

if 条件式:(;など)(改行)
　　実行文
else:(;など)
　  実行文

のようにif文などを書いていきます。

　if分では条件演算子を用いて数値や文字の比較をする式(条件式)を用いて、
その条件が真(True)だった時に配下に書いた実行文に移ります。
もし条件式が偽(False)だった時は次の文(elif,elseなど)に移ります。

最初の条件式のときにifを用い、他に分岐させたい項目がある場合はelifを用います。
ifやelifの条件から外れた物への指示はelse文の配下に記述します。

条件演算子には次のようなものがあります。

演算子	意味
>	左側の数が、右側の数よりも大きい
<	左側の数が、右側の数よりも小さい
==	両辺の数値の大きさが等しい
!=	等しくない
>=	左側の数が、右側の数よりも大きいか、または両辺が等しい
<=	左側の数が、右側の数よりも小さいか、または両辺が等しい

また条件を掛け合わせる時にはビット演算子を使うと便利です。

演算子	意味
&	両方とも成り立つ。
|	少なくとも、どちらか一方は正しい。

if_sequence.py

day = int(input("What's the date today ? "))
week_of_day = input("What day is today ? ")

if (day==13) & (week_of_day=="Friday") :
    print("Jason will come.")
else:
    print("Jason will not come.")

条件分岐は１行でも記述することができます。

score = int(input("score ?"))
print ("Good Score!" if score >= 70 else "Not Good Score...")

繰り返しをしてみよう

繰り返しにはfor文とwhile文の2種類があります。

for文

　for文は以下のような文法で成り立っています。
ここでデータの集合とは１〜９の数字の集合であったり、配列であったり、
ファイルの集合であるフォルダであったりします。

for 変数 in データの集合:
    実行文

実際にfor文を使ってみましょう。
今回はデータの集合に0~9の数字の集合を使うことで10回の操作を繰り返します。
0~9までの数字の集合には組み込み関数のrangeを用います。

count.py

for i in range(0,10):
    if i == 5:
        print("a")
        continue
    """if i == 9:
        break"""
    print(i)

i==5の時以外は0~9までの数字が出力されました。print(i)が動いてますね。
i==5の時だけaが出力されました。if文の中に入っていることがわかります。
if文の最後の行にcontinueと書いてあります。
coutinue文に当たると現在進行中の繰り返し処理一回分を終了させて
強制的につぎの繰り返し処理に入っていきます。
なのでi==5の処理時にはcontinue以降のprint(i)が実行されていません。

この他にもbreak文というのもあります。
これは繰り返し自体を終了させてしまいます。
i==9のコメントアウトを外して実行すると出力が8までになるはずです。

while文

while文の文法は以下のようになっています。

while 条件式:
    実行文

while文では条件式が真(True)である限り実行式が繰り返し実行されます。
for文と同じプログラムを書いていきましょう。

i = 0
while i < 10:
    if i == 5:
        print("a")
    elif i == 9:
        break
    print (i)
    i += 1

関数

　Pythonでは関数の定義も行うことができます。

f(x,y) = 2x + y

を出力するような関数を作ってみます。

def function(x=0,y=0):
    print("2x+y = %f"%(2*x+y))

x = int(input("x = "))
y = int(input("y = "))
function(x,y)

このように関数は

def 関数名(引数):
    実行文

のような文法で記述されます。引数は何個でも大丈夫です。
引数にはデフォルトの値を入れることができます。
関数を呼び出した時に引数の場所に値を入れていないとデフォルト値を入れた状態で実行されます。

返り値がある場合

return文を使うことで関数の返り値を設定できます。

def function(x,y):
    return 2 * x + y
x = int(input("x = "))
y = int(input("y = "))
print("2x+y = %f"%(f(x,y)))

課題

任意の値(n)までのフィボナッチ数列を出力する関数を作ってみましょう。

スコープ

　プログラムの変数、（メソッド、クラス名）には有効範囲があります。これをスコープと呼びます。
スコープは対象を宣言した場所によって変わります。
スコープは３種類に分けられます。

スコープ名	意味
ビルドインスコープ	Pythonプログラムの中に常に存在するスコープ
グローバルスコープ	作成するPythonファイルやプロジェクト内で有効
ローカルスコープ	関数内でのみ有効

スコープの大きさは、
ローカルスコープ⊂グローバルスコープ⊂ビルドインスコープ
になっています。スコープの外側からは、その名前を参照できません。

クラス

　クラスの構造は次のようになっています。
クラスとはデータ構造を作る仕組みで、クラスを使うと新しいデータ型を作ることができます。

class Personal_data():
    def __init__(self, name, age, weight, height): #コンストラクタの定義
        # インスタンス変数（属性の初期化）
        self.name = name
        self.age = age
        self.weight = weight
        self.height = height
        print("コンストラクタが呼ばれました")

    def hello(self):
        print("Hello, " + self.name)
    def BMI(self):
        print(self.weight/(self.height/100)**2)
    def __del__(self):                            #デストラクタの定義
        print('デストラクタが呼ばれました')

iwana = Personal_Data("Hiroki",23,65,170) # iwanaというインスタンスを生成しています。
iwana.hello()
iwana.BMI()
del iwana

　このようにクラスには変数や関数といった様々なものを入れることができます。
クラスの中に定義された関数をメソッド、クラスが内包しているデータを属性と呼びます。
クラスや関数は定義しただけでは何も行われません。
クラスや関数などを実際に使う状態にすることをインスタンスの作成と呼びます。
　コンストラクタはクラスのインスタンスが作成される際に一度だけ呼ばれるメソッドです。
今回ではiwana = Personal_Data("Hiroki",23,65,175)のときにParsonal_Dataクラスのインスタンスが作成され、そのオブジェクトを変数iwanaに代入しています。
オブジェクトとは属性とメソッドを持ったもののことを指します。
　コンストラクタとは別にクラスの特殊なメソッドとして、デストラクタがあります。デストラクタは呼び出された時にインスタンスを削除します。

課題

自分で好きなクラスを作ってみよう！

Pythonとは？

• Pythonの特徴

動的な型付け	型宣言をプログラム実行前にしなくて良い
ガベージコレクション	不要なメモリは自動的に解放
マルチパラダイム	オブジェクト指向、命令型、手続き型、関数型など複数の特性を併せ持つプログラミング言語
外部ライブラリが豊富	機械学習用の「TensorFlow」「Chainer」や Web開発用の「Django」「Flask」など様々
読みやすい	簡単な短いコードによって書けて読みやすい
多くのプラットフォームに対応	Windows,Mac,Linux,ios,Android...

C言語でのプログラム実行まで

たった１行の出力でもいろいろしなくちゃいけない...
＊もちろんC言語も高速であったり自由度が高かったりと優れているところがたくさんあります。

Hello.c

//ライブラリの読み込み
#include <stdio.h>
//main関数の定義
int main(int argc, char *args[]){
    printf("Hello, world!\n");
    return 0;
}

$  gcc -o hello hello.c　　　　　　　//コンパイル
$  ./hello　　　　　　　　　　　　　　 //実行
Hello, world!
$

Pythonでのプログラム実行まで

たったこれだけで動きます。
Pythonはインタプリタ言語といって
コンピュータがプログラミング言語を機械語に翻訳する作業を
プログラム実行時に同時作業で行ってくれます。
これにより、機械語への変換を一括で行うコンパイルという作業の指示をする必要がありません。

Hello.py

print("Hello, world!")

$ python Hello
Hello, world!

---

実際に書いてみよう！

Google Colaboratory
* ファイルをクリック

* ノートブックを新規作成をクリック

* 以下のコードを記入

print("Hello, world!")

• ボタンを押す

変数を使ってみよう

msg = "Hello, world!"
num = 100
ADMIN_EMAIL = "hello@mail.com"
squares = [1,4,9,16,25]
print(msg)
print(num)
print(squares)
print(ADMIN_EMAIL)


msg = "I'm Pretender"
num = num + 50
squares += [36, 49, 64, 81, 100] 
print(msg)
print(num)
print(squares)

Hello, world!
#numの出力
100
#squaresの出力
[1, 4, 9, 16, 25]
#ADMIN_MAILの出力
hello@mail.com
#更新後のmsgの出力
I'm Pretender
#更新後のmsgの出力
150
#更新後のsquaresの出力
[1, 4, 9, 16, 25, 36, 49, 64, 81, 100]

数値の四則演算など

x = int(input("x="))
print(x+3)
print(x-3)
print(x*3)
print(x/3)
msg += "GoodBye!"
print(msg)

文字列に値を埋め込む

name = "Iwana"
height = 173
#やり方１　文字列フォーマット指定
print("%sくんの身長は、%dcmです。"%(name,height)) 
#やり方２ {}で指定
print("{1}くんの身長は、{2}cmです。".format(name,height))

フォーマット指定子	型
%d	整数
%f	浮動小数
%s	文字列
%o	8進数
%x	16進数

フォーマット指定子のサイト

課題

名前とBMIを計算して表示するプログラムを書いてみよう。

条件分岐をしてみよう

if文

score=int(input("score?"))
if score >= 90:
    print("A+")
elif score >= 80:
    print("A")
elif score >= 70:
    print("B")
elif score >= 60:
    print("C")
else:
    print("You're failed.")

Pythonに限らずプログラミング言語では基本的に

if 条件式:
　　実行文
else:
　  実行文

のようにif文などを書いていきます。

　if分では条件演算子を用いて数値や文字の比較をする式(条件式)を用いて、
その条件が真(True)だった時に配下に書いた実行文に移ります。
もし条件式が偽(False)だった時は次の文(elif,elseなど)に移ります。

最初の条件式のときにifを用い、他に分岐させたい項目がある場合はelifを用います。
ifやelifの条件から外れた物への指示はelse文の配下に記述します。

条件演算子には次のようなものがあります。

演算子	意味
>	左側の数が、右側の数よりも大きい
<	左側の数が、右側の数よりも小さい
==	両辺の数値の大きさが等しい
!=	等しくない
>=	左側の数が、右側の数よりも大きいか、または両辺が等しい
<=	左側の数が、右側の数よりも小さいか、または両辺が等しい

また条件を掛け合わせる時にはブール演算子を使うと便利です。

演算子	意味
and	両方とも成り立つ。
or	少なくとも、どちらか一方は正しい。

if_sequence.py

day = int(input("What's the date today ? "))
week_of_day = input("What day is today ? ")

if (day==13)and(week_of_day=="Friday") :
    print("Jason will come.")
else:
    print("Jason will not come.")

条件分岐は１行でも記述することができます。

score = int(input("score ?"))
print ("Good Score!" if score >= 70 else "Not Good Score...")

繰り返しをしてみよう

繰り返しにはfor文とwhile文の2種類があります。

for文

　for文は以下のような文法で成り立っています。
ここでデータの集合とは１〜９の数字の集合であったり、配列であったり、
ファイルの集合であるフォルダであったりします。

for 変数 in データの集合:
    実行文

実際にfor文を使ってみましょう。
今回はデータの集合に0~9の数字の集合を使うことで10回の操作を繰り返します。
0~9までの数字の集合には組み込み関数のrangeを用います。

count.py

for i in range(0,10):
    if i == 5:
        print("a")
        continue
    """if i == 9:
        break"""
    print(i)

i==5の時以外は0~9までの数字が出力されました。print(i)が動いてますね。
i==5の時だけaが出力されました。if文の中に入っていることがわかります。
if文の最後の行にcontinueと書いてあります。
coutinue文に当たると現在進行中の繰り返し処理一回分を終了させて
強制的につぎの繰り返し処理に入っていきます。
なのでi==5の処理時にはcontinue以降のprint(i)が実行されていません。

この他にもbreak文というのもあります。
これは繰り返し自体を終了させてしまいます。
i==9のコメントアウトを外して実行すると出力が8までになるはずです。

while文

while文の文法は以下のようになっています。

while 条件式:
    実行文

while文では条件式が真(True)である限り実行式が繰り返し実行されます。
for文と同じプログラムを書いていきましょう。

i = 0
while i < 10:
    if i == 5:
        print("a")
    elif i == 9:
        break
    print (i)
    i += 1

関数

　Pythonでは関数の定義も行うことができます。

f(x,y) = 2x + y

を出力するような関数を作ってみます。

def function(x=0,y=0):
    print("2x+y=%f"%(2*x+y))

x = int(input("x="))
y = int(input("y="))
function(x,y)

このように関数は

def 関数名(引数):
    実行文

のような文法で記述されます。引数は何個でも大丈夫です。
引数にはデフォルトの値を入れることができます。
関数を呼び出した時に引数の場所に値を入れていないとデフォルト値を入れた状態で実行されます。

返り値がある場合

return文を使うことで関数の返り値を設定できます。

def function(x,y):
    return 2*x + y
x = int(input("x="))
y = int(input("y="))
print("2x+y=%f"%(f(x,y)))

課題

任意の値(n)までのフィボナッチ数列を出力する関数を作ってみましょう。

スコープ

　プログラムの変数、（メソッド、クラス名）には有効範囲があります。これをスコープと呼びます。
スコープは対象を宣言した場所によって変わります。
スコープは３種類に分けられます。

スコープ名	意味
ビルドインスコープ	Pythonプログラムの中に常に存在するスコープ
グローバルスコープ	作成するPythonファイルやプロジェクト内で有効
ローカルスコープ	関数内でのみ有効

スコープの大きさは、
ローカルスコープ⊂グローバルスコープ⊂ビルドインスコープ
になっています。スコープの外側からは、その名前を参照できません。

クラス

　クラスの構造は次のようになっています。
クラスとはデータ構造を作る仕組みで、クラスを使うと新しいデータ型を作ることができます。

class Personal_data():
    #コンストラクタの定義
    def __init__(self, name, age, weight, height): 
        # インスタンス変数（属性の初期化）
        self.name = name
        self.age = age
        self.weight = weight
        self.height = height
        print("コンストラクタが呼ばれました")

    def hello(self):
        print("Hello, " + self.name)
    def BMI(self):
        print(self.weight/(self.height/100)**2)

    #デストラクタの定義
    def __del__(self):                            
        print('デストラクタが呼ばれました')
# iwanaというインスタンスを生成しています。
iwana = Personal_Data("Hiroki",23,65,170) 
iwana.hello()
iwana.BMI()
del iwana

　このようにクラスには変数や関数といった様々なものを入れることができます。
クラスの中に定義された関数をメソッド、クラスが内包しているデータを属性と呼びます。
クラスや関数は定義しただけでは何も行われません。
クラスや関数などを実際に使う状態にすることをインスタンスの作成と呼びます。
　コンストラクタはクラスのインスタンスが作成される際に一度だけ呼ばれるメソッドです。
今回ではiwana = Personal_Data("Hiroki",23,65,175)のときにParsonal_Dataクラスのインスタンスが作成され、そのオブジェクトを変数iwanaに代入しています。
オブジェクトとは属性とメソッドを持ったもののことを指します。
　コンストラクタとは別にクラスの特殊なメソッドとして、デストラクタがあります。デストラクタは呼び出された時にインスタンスを削除します。

課題

自分で好きなクラスを作ってみよう！

はじめに

　新型コロナウイルス感染症(COVID-19)の感染の中心地はもはや中国本土から欧米に移行しつつあり、EU域内でも自由な国境の移動を制限したり、外出禁止令が発令されるなどの緊急事態になっています。日本において、新型コロナウイルスの脅威が認識されはじめたのは、2020年2月3日のダイヤモンド・プリンセス号の横浜港入港、2月25日の政府の基本方針表明、2月28日の北海道における緊急事態宣言、などが契機ですが、具体的な対策として2月26日から大規模なイベントの自粛要請、3月2日からの全国の小中高校の臨時休校などがなされて来ました。
　それから約半月が経過しますが、その間、日経平均株価が22000円台から17000円台まで急落するなど、先行きに暗雲が立ち込めています。新型コロナウイルス感染症の影響は健康のみならず、社会経済に深刻な打撃を与えており、いつ収束するのかが極めて重大な関心を集めています。
　そこで、本記事では、１月末から３月頭にかけて、日本における新型コロナウイルスの感染力の指標である、基本再生産数がどのように推移してきたか、計算で求める一つの方法を提示したいと思います。
　計算式やコードは少し長いので、計算結果を先に見て頂いた方が良いかもしれません。

基本再生産数とは？

　基本再生産数は、実は非常に奥が深いようで、京大数理解析研究所講究録や日本数理生物学会レター等を見ると、世代間の変遷を考慮するなど単純ではないようですが、敢えて、一言でいうと、

• 人口集団の典型的な１人の感染者が、その全感染期間において再生産する２次感染者数の期待値

と定義されるようです。この数値をR0とすると、感染が拡大する条件は$R_0>1$、感染が収束する条件は$R_0 < 1$となります。

　イギリスのボリス・ジョンソン首相が、3月16日の会見で、イギリス国民の約60%が集団免疫を獲得するまでピークをコントロールする趣旨の説明をしたようですが（ニューズウィークの記事）、これは、$R_0 = 2.6$と仮定して、ワクチンや自己免疫などでの集団免疫率$H$が、

(1-H)R_0 < 1

を満たす条件から$H > 1-1/R_0=0.615$と想定されたようです。つまり、集団免疫によって、$R_0$の条件が緩和されるのですね。しかし、ワクチンの開発には臨床試験など１年以上かかる見込みですし、敢えて過剰に感染を抑えないことで集団免疫の獲得をしようという大胆な戦略です。

基本再生産数の計算モデル

　ここで、基本再生産数を計算するために、上記の定義とSEIRモデルをベースとした、簡易モデルを考えたいと思います。次の図を見てください。

上段は日付、下段の棒はある日の感染者集団を表しています。各記号の意味は、

• E:感染症が潜伏期間中の者(Exposed)
• I:発症者(Infectious)
• P:検査陽性者(Positive)
• lp:潜伏期間
• ip:感染期間

を表しています。つまり、ある日tにおけるIからEの生産は、t以前に感染し発症した状態にあるIの集団からの感染であり、その日tの再生産数は$R_0(t)/ip$である、というモデルです。
　これを式で表すと、

\frac{1}{ip} R_0(t) \times \sum_{s=t+1}^{t+ip} P(s) = P(t+lp+ip), \\
\therefore R_0(t) = \frac{ ip \times P(t+lp+ip)} {\sum_{s=t+1}^{t+ip} P(s)}

となります。以下、この式を使って、公開されている検査陽性者のデータから、基本再生産数$R_0(t)$の時間推移を計算してみましょう。

元になるデータ

　本記事では、都道府県別新型コロナウイルス感染者数マップ（ジャッグジャパン株式会社提供）において、公開されているcsvデータを使わせて頂きました。都道府県のデータが１つに集約されており、大変使いやすい形になっています。

Pythonで計算してみる

では、Pythonを使って、都道府県別の基本再生産数R0(t)を計算してみましょう。

前提条件

前提条件として、以下を用います。

• lp:潜伏期間 = 5[day]
• ip:感染期間 = 8[day]
• 上記の図のPは確定診断日とする。

なお、lpとipには諸説ありますが、一つの目安として設定しています。また、計算式の都合上、データの最新の日付から、lp+ip=13[day]以前の$R_0(t)$しか計算できません。グラフ上では、13日前～最新日の期間の値が0になっていても意味のない値ですので、ご注意ください。

コード

まず、ライブラリをインポートします。

# coding: utf-8

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import datetime
import locale

都道府県別新型コロナウイルス感染者数マップ（ジャッグジャパン株式会社提供）で提供されているcsvから、各都道府県別のデータを読み込み、データフレームに変換する関数です。

def readCsvOfJapanPref(pref : None):
    # 下記URLよりダウンロード
    # https://jag-japan.com/covid19map-readme/
    fcsv = u'COVID-19.csv'
    df = pd.read_csv(fcsv, header=0, encoding='utf8', parse_dates=[u'確定日YYYYMMDD'])
    # 確定日, 受診都道府県のみ抽出
    df1 = df.loc[:,[u'確定日YYYYMMDD',u'受診都道府県']]
    df1.columns = ['date','pref']
    # 受診都道府県で抽出
    if pref is not None: # Noneであれば日本全体
        df1 = df1[df1.pref == pref]
    df1 = df1.loc[:,'date']
    # 確定日でカウント
    df2 = pd.DataFrame( df1.value_counts() ).reset_index()
    df2.columns = ['date','P']
    df2 = df2.sort_values('date')
    return df2

計算用のデータフレームを定義する関数です。列のPpreは計算式の分母用、Patは計算式の分子用です。2020年1月24日を開始日として、days分だけ計算の枠を作ります。

def makeCalcFrame(days):
    t_1 = pd.Timestamp(2020,1,24) # 計算開始日
    td = pd.Timedelta('1 days')
    #
    npd = [[t_1 + td * i, 0, 0, 0 ] for i in range(0,days)]
    df1 = pd.DataFrame(npd)
    df1.columns = ['date', 'Ppre','Pat', 'R0']
    #
    return df1

都道府県データと計算フレームを結合する関数です。

def mergeCalcFrame(df1, df2):
    return pd.merge(df1, df2, on='date', how='left').fillna(0)

計算式に従って、$R_0$を計算する関数です。計算を簡単にするため、インデックスの範囲外にアクセスしたらNaNを返すラムダ式を作りました。

def calcR0(df, keys):
    lp = keys['lp']
    ip = keys['ip']
    nrow = len(df)
    getP = lambda s: df.loc[s, 'P'] if s < nrow else np.NaN
    for t in range(nrow):
        df.loc[t, 'Ppre'] = sum([ getP(s) for s in range(t+1, t + ip + 1)])
        df.loc[t, 'Pat' ] = getP(t + lp + ip)
        if df.loc[t, 'Ppre'] > 0:
            df.loc[t, 'R0'  ] = ip * df.loc[t, 'Pat'] / df.loc[t, 'Ppre']
        else:
            df.loc[t, 'R0'  ] = np.NaN
    return df

結果をグラフ表示する関数です。

def showResult(df, title):
    # R0=1 : 収束のためのターゲット
    ptgt = pd.DataFrame([[df.iloc[0,0],1],[df.iloc[len(df)-1,0],1]])
    ptgt.columns = ['date','target']
    # show R0
    plt.rcParams["font.size"] = 12
    ax = df.plot(title=title,x='date',y='R0', figsize=(10,7))
    ptgt.plot(x='date',y='target',style='r--',ax=ax)
    ax.grid(True)
    ax.set_ylim(0,)
    plt.show()

各都道府県に対する計算を行う部分です。

def R0inJapanPref(pref, label):
    keys = {'lp':5, 'ip':8 }
    df1 = makeCalcFrame(60) # 60 days
    df2 = readCsvOfJapanPref(pref)
    df = mergeCalcFrame(df1, df2)
    df = calcR0(df, keys)
    showResult(df, 'COVID-19 R0 ({})'.format(label))
    return df

preflist = [[None, 'Japan'], [u'東京都', 'Tokyo'],\
            [u'大阪府', 'Osaka'],  [u'愛知県', 'Aichi'],\
            [u'北海道', 'Hokkaido']]
dflist = [[R0inJapanPref(pref, label), label] for pref, label in preflist]

上の計算結果をまとめて表示する部分です。

def showResult2(ax, df, label):
    # show R0
    plt.rcParams["font.size"] = 12
    df1 = df.rename(columns={'R0':label})
    df1.plot(x='date',y=label, ax=ax)

# R0=1
ptgt = pd.DataFrame([[df_tko.iloc[0,0],1],[df_tko.iloc[len(df_tko)-1,0],1]])
ptgt.columns = ['date','target']
ax = ptgt.plot(title='COVID-19 R0', x='date',y='target',style='r--', figsize=(10,8))
#
for df, label in dflist:
    showResult2(ax, df, label)
#
ax.grid(True)
ax.set_ylim(0,10)
plt.show()

計算結果

　それでは、いよいよ計算結果を見てみましょう。

日本全体の基本再生産数の推移

• 2月10日までは最大10付近までの大きな流入があるようです。これは、恐らく、春節（1月24日～1月30日）に来日した中国旅行者の影響ではないでしょうか。
• その後は、一旦3付近まで上昇するものの、1を辛うじて超える付近（平均して1.5付近）で推移しているようです。

東京都の基本再生産数の推移

• 2月6日ごろまでは最大16の高い値でしたが、2月7日以降はおおむね落ち着いているようです。
• ただし、2月18日～2月23日や3月1日など、時折3を超える山があり、クラスターの発生や、海外からの流入があった可能性を示唆しています。

大阪府の基本再生産数の推移

• 2月17日～2月24日に最大26の大きなピークが見られます。これは恐らくクラスターの発生でしょう。
• 実際、高槻市のページによると、大阪市内の４つのライブハウスで2月15日から2月24日にかけて行われたライブの参加者から、感染クラスターの発生が確認されています。

愛知県の基本再生産数の推移

• 2月12日までの、恐らく海外渡航者による感染が見られます。
• その後、2月18日から2月26日にかけておおむね3付近で推移する山が見られます。これは、報道されているように、ジムや福祉施設で起きたクラスターと思われます。
• クラスターの山と山が隣接してしまっているので、クラスター連鎖が心配されます。

北海道の基本再生産数の推移

• 2月8日から2月10日にかけて最大80もの鋭いピークが見られます。これは恐らく、さっぽろ雪祭り（2月4日～2月11日）で発生したクラスターではないかと思われます。
• 2月16日以降は、概ね1を下回る水準で推移しており、収束に向かっているものと思われます。

全体を比較

• 2月18日までに主に海外からの流入を主とした感染、3月1日までに国内でのクラスターの散発的な発生があったものと思われます。
• 3月1日以降は、全般的に落ち着いてきており、まだ1を超える水準であるものの、ギリギリ収束に向かっていると言えるのではないでしょうか。

考察

以上から、日本の各都道府県の感染者データに基づく、基本再生産数R0の推定に関して、次の傾向がシミュレーションから導けます。

• 2月中旬までの海外からの流入が最もR0が高く、その後国内でのクラスターと思われる影響で散発的に跳ね上がる傾向が見られる。
• それ以外の期間では、特に北海道では1を下回る期間をキープ出来ており、概ね良好な状態。
• 計算結果のクラスター疑惑日と、報道されているクラスター発生日（大阪のライブハウス、さっぽろ雪まつり）がよく一致しており、計算モデルからクラスター発生を当てられる可能性が高い。（もっとも、積極的疫学調査による発見の効果もあると思われます。）

さらに言えば・・・

• シミュレーション結果からは、基本再生産数$R_0$は一定ではなく、かなり変動しています。
• 基本再生産数$R_0$を上昇させる主要因は、海外からの流入と、クラスター発生であることが数値的にも裏付けられます。これらを抑制する対策は、やはり海外渡航制限と、クラスターの要因になりうるイベントや施設の自粛であると思われます。
• 全体的な傾向では、3月1日以降は収束傾向にあり、日本全体でR0が1に近づいているため、現状の対策が功を奏しているようにも見えます。
• $R_0$は陽性者数と陽性者数の比率なので、検査されていない隠れ感染者が検査陽性者の一定倍いたとしても影響は軽微でしょう。

参考リンク

下記のページを参考にさせて頂きました。
都道府県別新型コロナウイルス感染者数マップ（ジャッグジャパン株式会社提供）
京大数理解析研究所講究録
日本数理生物学会レター
ニューズウィークの記事
SEIRモデル
高槻市のページ

直接引用してませんが、大変有用なリンクを追加しておきます。
COVID-19 reports
全国クラスターマップ

概要

自作の Middleware を試験したいけど、別の Middleware に依存しているケースの対処法です。

どういうケース？

一般的に Middleware をテストする場合、Middleware のインスタンスを作成し、process_request メソッドを呼ぶのがメジャーなようです。
参照:Unit Testing Django Middleware - Adam Donaghy - Medium

しかし、別の Middleware で何かしらの処理を行っていることを前提にしている Middleware について、上記の方法では対応できません。
例えば、認証情報に係る処理を行う Middleware の場合、セッション関連のミドルウェアなどを経由していないと request.user 参照時に user attributeが存在せず例外が発生することなどが考えられます。

対処

• self.client.get() で実行します (View を試験するのと同じように)
• Middleware が下の層を呼び出す get_response() を適当なMockに差し替えます
  
  • ※クラスインスタンスが持つ get_response() は __init__ で注入するものになっているので、__init__ を適当なメソッドにすり替える必要があります。

サンプル

Middleware

class SampleMiddleware(object):
    def __init__(self, get_response):
        self.get_response = get_response

    def __call__(self, request):
        # いろいろやる

        return self.get_response(request)

テスト

def fake_init(self, _):
    def get_response(request):
        return HttpResponse()

    self.get_response = get_response


class SampleMiddlewareTest(TestCase):
    def setUp(self):
        self.mock_init = mock.patch.object(
            AuthMiddleware, "__init__", fake_init
        ).start()

        self.addCleanup(mock.patch.stopall)

    def test_いろいろテスト(self):
        # リクエストする
        self.client.get("/")

        # assertする

参照

• Mocking complicated __init__ in Python - George Shuklin - Medium
  • __init__ を lambda で何もしない関数を注入しているサンプルを参考にしました。