投資のためのデータサイエンス

個人の投資活動に役立つデータ分析にまつわる話題を綴ります。

データサイエンスのおすすめオンライン記事(10月20日付)

先週「ワーケーションツアー」なるものに参加した。様々な所から集まった参加者が観光を楽しみつつ、ある時間帯にはPCに向かって仕事をするのは異様な感じもしたが、休暇を楽しみつつ会社ネットワークにつないで「置いてきぼり」にならないようにするにはいい方法だとも思った。他の参加者との交流もできて有意義だった。

今週も7月~8月ごろに新たに掲載された記事をピックアップして紹介する。

Pythonを用いて無料で貴方の画像をラベル付けする簡単な方法
towardsdatascience.com
画層認識で重要になる画像のラベル付けに関するチュートリアル記事。

貴方のノートブックをJupyterLabテンプレートに同期させる
towardsdatascience.com
Jupyter Notebook上でコード開発をすると前に作ったコードのコピーペーストの連続になりがちだが、ここではよりスマートな方法が紹介されている。

データサイエンスをリードする10人の女性
www.siliconrepublic.com
2020年10月の記事。10人の女性データサイエンティストの紹介。

重要なデータサイエンスのソフトスキル:成長マインドセット
benjaminobi.medium.com
データサイエンティストとして必要なソフトスキル(定性的・人間的スキル)のうち、「成長マインドセットを持ち続ける」ための具体的な方策が述べられている。

データサイエンスの将来トレンド2022
www.datasciencecentral.com
データサイエンスの今後のトレンドについて、主に政策的視点から述べられている。「AIの民主化」「トッププログラミング言語は依然Python」など。

データサイエンスのおすすめオンライン記事(9月28日付)

事実上の次期首相を選ぶ選挙が翌日に迫ってきた。次のリーダーには是非長期的な視点で日本の経済的衰退を食い止める政策を打ってもらいたい。特にエネルギー政策は国の根幹に関わる重要課題であり、次期首相の政治家としての手腕に期待したい。
前回の投稿からかなり間が空いてしまい、記事のストックが積み上がってしまっている。今回紹介する記事も6月~7月ごろの記事が中心である。

貴社のCRMデータのよくある問題と解決法
www.datasciencecentral.com
Customer Relationship Management(顧客関係管理)データを扱う場合によく遭遇する問題とその解決策をまとめている。

データサイエンティストとして貴方の再利用可能なPythonコードを管理する
www.kdnuggets.com
プログラムコードの管理と再利用は仕事でコードを書く人なら誰でも直面する問題であろう。ここでは初心者向けに個人のコード再利用の方法が述べられている。

「バリュー」の定義:AI成功のカギ
www.datasciencecentral.com
おなじみ"Dean of Big Data"ことBill Schmarzo氏のエッセイ。「データサイエンティストのように考える」フレームワークの改訂版と「AIプロジェクトがもたらすバリューを定義する重要性を述べている。

Facebookの"Prophet"は時系列予測の救世主か、それともただのわんぱく小僧か
www.datasciencecentral.com
Facebookが開発し、最近注目を集めている時系列予測ツールProphetについての解説記事。

最新のAI研究論文を手早く分析する必須ツール
www.datasciencecentral.com
AIの最新論文をリサーチするのは骨のおれる作業だが、ここではData Science Centralによく寄稿している著者が利用している"labml.ai"というプラットフォーム上の論文サーチ機能を紹介している。

Pythonによる分散分析のデモ

ある変量の変動の要因を分析する場合には、母集団のサブグループ別に分布や平均値を比較して、グループ間で差があるかどうかをみる、ということはよく行われます。これをより系統的に行う手法が分散分析です。分散分析は、原理的にはダミー変数を用いた回帰分析と同じですが、例えば「性別」や「年齢階級」といった要因単位で指標値に有意に影響を及ぼしているかを見るような場合に便利です。ここでは都道府県別の各種指標値を例題として、「(地理的に隣接した)地方」や「大都市があるか」といった要因で異なるかどうかを見ることにします。
まず、必要なライブラリをインポートします。

# ライブラリのインポート
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
sns.set(font=['IPAexGothic'])
import japanize_matplotlib # matplotlib日本語
from google.colab import files # 画像ファイルダウンロード用

ここでは、日本の47都道府県をクロスセクションデータとして用いることにします。まず、エクセルファイルからPandasデータフレームに都道府県の属性データを読込みます。

# ファイルからデータを読み込む
df01 = pd.read_excel('/content/drive/My Drive/Colab Notebooks/データ解析講座一日目演習.xlsx',sheet_name='都道府県')
df01.head()

f:id:nicjps230:20210902182435p:plain
上記のデータ列のうち、「都市化」は政令指定都市がある道府県から北海道を除いて東京都と香川県(四国経済の拠点が多く存在)を加えたものを「都市県」、その他を「郊外県」としました。また、「メディア」は以下の記事にあるメディア接触の特徴に関する分類です。 https://news.livedoor.com/article/detail/9064551/
同様に、公的ウェブサイトに掲載されている、47都道府県の家計指標データを読込みます。

# ファイルからデータを読み込む
df02 = pd.read_excel('/content/drive/My Drive/Colab Notebooks/データ解析講座一日目演習.xlsx',sheet_name='家計指標')
df02.head()

f:id:nicjps230:20210902182854p:plain
そして、読み込んだ2つのデータフレームを、県名をキーとして横に連結します。

# 県名をキーにして横連結する
df1 = pd.merge(df01, df02, left_on='県名', right_on='県名')
df1.drop(columns=['コード_x','コード_y'],inplace=True)

一元配置分散分析
ここでは都道府県の食料費割合の変動に及ぼす「地域」等の属性要因の影響をみるために、一元配置分散分析を実行します。statsmodelsライブラリを用います。

# 一元配置分散分析
import statsmodels.formula.api as smf
import statsmodels.api as sm
oneway_anova = smf.ols(formula='食料費割合~地域+都市化+メディア',data=df1).fit()
sm.stats.anova_lm(oneway_anova,typ=2).round(2)

これにより、以下のように各要因についてのF検定統計量が出力されます。
f:id:nicjps230:20210902183641p:plain
尚、交互作用項を加える場合はモデル式に"地域*都市化"のように*でつないだ項を追加します。
有意な主効果の様子を見る
ここではseabornライブラリのcatplot()を用いて、各水準毎の指標値の平均を棒グラフで比較します。

# 有意差のある効果についてその中身をみる
import matplotlib.pyplot as plt
import seaborn as sns
p = sns.catplot(x='食料費割合',y='都市化',data=df1,kind='bar',height=4,aspect=2)
p.set(title='食料費割合の都市化による差')

以下のような形で平均の違いを比較することができます。
f:id:nicjps230:20210902184340p:plain
交互作用の様子を見る
交互作用(ある要因による指標への影響が、他の要因の水準によって異なるかどうか)の様子を見るためには、statsmodelsライブラリのinteraction_plotを用います。

from statsmodels.graphics.factorplots import interaction_plot
fig=interaction_plot(df1.地域,df1.都市化,df1.食料費割合,colors=['red','blue'],ms=10)
plt.title('地域別食料費割合と都市化')
plt.legend(bbox_to_anchor=(1.05,1.0),loc='upper left')

これにより交互作用の様子がひと目でわかるようになります。近畿地方において、都市部と郊外部で食料費割合の差が大きいことがわかります。
f:id:nicjps230:20210902184835p:plain

データサイエンスのおすすめオンライン記事(8月9日付)

前回のこのシリーズの投稿からバタバタしている間にオリンピックが終わってしまった。これは最近始まったことではないが、バレーボールではコーチがタブレット端末片手にデータ分析結果に基づいてリアルタイムに戦略を指南していたし、サッカーでも日本代表チームには優れた分析スタッフがいて、自チームや対戦する相手チームの直近の試合のデータを徹底的に分析しているという話をしていた。スポーツのデータ活用は予算と人材さえあればどこでも始められるようになってきている。
今回は「機械学習」と「スクレイピング」の記事を中心に紹介する。

Pythonコードを高速化する方法
www.kdnuggets.com
インタープリター言語であるPythonを少しでも早く走らせるプログラミング上のコツをいくつか挙げている。

ウェブページのスクレイピングに回帰ベースの機械学習モデルを用いる
www.datasciencecentral.com
ユーザエージェントとブロックされる可能性の関連をモデル化してスクレイピングの効率化を図る。

225個のPythonによる機械学習プロジェクト
medium.com
Python言語を用いた機械学習プロジェクト(問題・データ・コードが掲載されており自分で試して学べる)が225個も紹介さている良記事。初級者向けプロジェクトは最初の方でまとめて紹介されている。

Amazonの商品カタログデータをスクレイピングする方法
www.datasciencecentral.com
Pythonのscrapyライブラリを用いてAmazonの商品カタログデータをスクレイピングする具体的な方法の解説。

顔認識:最新の方法と最良のツール
www.datasciencecentral.com
「知識ベース」「テンプレート照合」「特徴量ベース」「外見ベース」といった顔認識の最新技術の解説と有料/無料で利用できる顔認識ツールの紹介。

PythonによるECカタログデータの分析小技集(下)

前回に引き続き、「楽天ブックス書籍検索API」を利用してダウンロードした書籍情報データを用いて、ECサイトのカタログデータを処理・分析するためのpythonの小技をまとめました。
まず前回と同様に必要なライブラリを読み込み、データを保存先からロードします。

# ライブラリのインポート
import numpy as np
import pandas as pd
import re
import matplotlib.pyplot as plt
import seaborn as sns
sns.set(font=['IPAexGothic'])
# parquetファイルからロード
df01 = pd.read_parquet('rakuten_books.parquet')
# 列を絞り込む
focus_cols=['title','subTitle','seriesName','contents','author','publisherName','size',\
	'itemCaption','salesDate','itemPrice','reviewCount','reviewAverage']
df1 = df01.loc[:,focus_cols]

定型テキストの出現頻度を求める
テキスト情報の列の中でも、出現する文字列のパターンが限定されている場合があります。いわゆる「カテゴリカルデータ」に相当するものです。その出現頻度を知りたいことはよくあります。これにはcollectionsライブラリのCounter()メソッドを用います。ここでは列'publisherName'(出版社)の頻度分布を求め、上位10社を出力します。

# 頻度分布を求める列名を指定する
col_name = 'publisherName'
import collections
c_count = collections.Counter(df1[col_name])
# カウント上位のリストを得る
p_names, p_counts = zip(*c_count.most_common(10))
print(p_names)

尚、出力結果は以下のようになりました。
('講談社', 'KADOKAWA', '集英社', '小学館', '宝島社', 'スクウェア・エニックス', '白泉社', '秋田書店', '幻冬舎', '朝日新聞出版')

正規表現による検索
決まった文字列ではなく、ある規則性を持った文字列を検索したい場合には「正規表現」による検索を用います。ここでは、コミックセットのように複数の書籍が一括して販売されているケースを検索し、「○冊セット」の○にあたる数字を取り出すことを考えます。ここではpythonライブラリ're'(regular expressionを略したもの)のsearchメソッド(と内容を取り出すためのgroupメソッド)を一件につき2回用いています。1回目の適用で'○冊セット'の文字列をとりだし、2回目の適用でその中の○にあたる数字の部分を取り出して、int()で数値化しています。

def func_bulk(s):
    """タイトル文字列sの中の冊数を表すパターンを検索して冊数の数値を返す"""
    nitem = 0
    blkstr1 = re.search('[0-9]+巻セット',s)
    if not blkstr1 is None:
        nitem = int(re.search('[0-9]+',blkstr1.group()).group())
        return nitem
    blkstr2 = re.search('[0-9]+冊セット',s)
    if not blkstr2 is None:
        nitem = int(re.search('[0-9]+',blkstr2.group()).group())
        return nitem
    blkstr3 = re.search('全[0-9]+巻',s)
    if not blkstr3 is None:
        nitem = int(re.search('[0-9]+',blkstr3.group()).group())
        return nitem
    blkstr4 = re.search('全[0-9]+冊',s)
    if not blkstr4 is None:
        nitem = int(re.search('[0-9]+',blkstr4.group()).group())
        return nitem
    return nitem
# セット書籍の冊数を取り出す
df1 = df1.assign(n_bulk = df1['title'].apply(func_bulk))
df81 = df1[df1['n_bulk'] > 0]
dfi.export(df81,'kensaku_bulk.png')

結果は以下のようになりました。
f:id:nicjps230:20210723155132p:plain

データフレーム内の複数の列の演算により新たな列を作成する
これまでにもデータフレーム内のある一つの列から、ある演算によって新しい列を作り出す例が出てきました。演算が簡単なものならばその列を作り出す文の中に式を埋め込めばいいのですが、式が複雑な場合は別途関数を定義して、元の列に関してapplyメソッドを使う、という形でした。それでは、データフレーム内の複数の列を使って新しい列を作り出す場合はどうするのでしょうか。その場合は、データフレームそのものについてappyメソッドを適用します。ここでは列'reviewCount'と列'reviewAverage'を用いて、新しい列'review_pi'を作成しています。

# 新しい列を作成するための関数
def func_review_pi(df):
  return np.log(df.loc['reviewCount']+1) + df.loc['reviewAverage']
# データフレーム内の複数の列の演算により新たな列を作成する
df51 = df1.assign(review_pi = df1.apply(func_review_pi, axis=1))

新しい列を作成するための関数に他の引数がある場合
これまでの例では、データフレーム内の新たな列を作るための関数には「元の列」や「元のデータフレーム」以外の引数はありませんでした。今度は、その関数が別途他の引数を必要とする場合の記述を見てみます。 ここでの例は、書籍のタイトルにある指定したキーワードが含まれていたら、あるジャンルの書籍であると判定するための関数の作成です。一つのジャンルについてキーワードは複数あるものとし、キーワードは2次元のリストで与えられるとします。

# ジャンル別のキーワードのリストの定義
# 0:経済, 1:歴史, 2:科学, 3:医学, 4:旅行, 5:辞典, 6:人生, 7:食事, 8:情報, 9:入試, 10:漫画
kws = [['経済','財務','お金','家計','会計','金融'],
       ['史','古代'],
       ['科学','サイエンス','宇宙'],
       ['医','癌','内科','外科','歯科','手術','疾患'],
       ['旅','紀行'],
       ['辞典','広辞苑'],
       ['人生','生き方','生きかた','ライフ'],
       ['料理','食事','グルメ','レシピ'],
       ['Excel','Word','Python','iPhone','iPad','プログラミング','アプリ','HTML','ソフトウェア','情報技術','Adobe'],
       ['入試','中学受験','高校受験','大学受験'],
       ['鬼滅の刃','リベンジャーズ','るろうに剣心','スライムだった件','青の祓魔師','呪術廻戦','進撃の巨人']]
# キーワード判定関数の定義
def func_genre(s, ww, k):
  """文字列sに2次元リストwwの第kジャンルのキーワードが
  含まれていれば1, 含まれていなければ0を返す"""
  q = 0
  for j in range(len(ww[k])):
    if ww[k][j] in s:
      q = 1
  return q

演算の元になる列以外の引数がある場合、その引数はapply()メソッドの中で引数名=値の形で指定します。

# タイトルに経済分野のキーワードを含むかどうかのダミー変数を生成する
df1= df1.assign(G経済=df1['title'].apply(func_genre, ww=kws, k=0))

PythonによるECカタログデータの分析小技集(上)

一部のECサイトでは、販売商品のカタログデータをAPIにより取得することができます。このようなデータには、商品名、商品の補足説明、価格などの情報が含まれています。ここでは、このようなECサイトのカタログデータを処理・分析するためのpythonの小技をまとめました。ここで用いるデータは「楽天ブックス書籍検索API」を利用してダウンロードした書籍情報データです。
まず必要なライブラリを読み込み、データを保存先からロードします。

# ライブラリのインポート
import numpy as np
import pandas as pd
import re
import matplotlib.pyplot as plt
import seaborn as sns
sns.set(font=['IPAexGothic'])
# parquetファイルからロード
df01 = pd.read_parquet('rakuten_books.parquet')
# 列を絞り込む
focus_cols=['title','subTitle','seriesName','contents','author','publisherName','size',\
	'itemCaption','salesDate','itemPrice','reviewCount','reviewAverage']
df1 = df01.loc[:,focus_cols]

テキスト項目の文字列検索
まず、テキスト情報からあるキーワードを含む行を抽出するケースを考えます。pandasのSeries.str.contains()を用います。これにより、指定した文字列を含む行のみを取り出すことができます。またこのコードでは、Jupyter Notebookの出力イメージをそのまま画像ファイルに変換するdataframe_imageライブラリを用いています。

import dataframe_image as dfi
# テキスト項目で特定の文字列を検索する
col_name = 'title'
ch_kw = '資産'
df201 = df1[df1[col_name].str.contains(ch_kw)]
dfi.export(df201,'kensaku_sisan.png')

以下のような検索結果が得られます。
f:id:nicjps230:20210709161416p:plain
次に、2つのキーワードについて、その両方を含む行を取り出したい(and条件)場合は、一つのキーワードで絞りこんだデータフレームに対してもう一つのキーワードでさらに絞りこみをかけることでできます。【注】一つの命令文にブール演算を埋め込んで検索することは、下に記載したquery()メソッドではできますが、それ以外ではエラーが出るためここでは避けています。

# テキスト項目で特定の複数の文字列をand検索する
col_name = 'itemCaption'
ch_kw1 = '成功'; ch_kw2 = '金'
df202 = df1[df1[col_name].str.contains(ch_kw1)]
df203 = df202[df202[col_name].str.contains(ch_kw2)]
print('該当するデータ件数 = ',len(df203))

一方、2つのキーワードのどちらかが含まれる行を抽出したい(or条件)場合は、元のデータフレームを個別のキーワードで絞りこんだ結果を縦結合することでできます。これは検索対象の列が異なる場合にも応用できます。どちらかの検索結果がゼロ行でも特に問題はありません。

# テキスト項目で特定の複数の文字列をor検索する
col_name = 'title'
ch_kw1 = '資産'; ch_kw2 = '財産'
df204 = df1[df1[col_name].str.contains(ch_kw1)]
df205 = df1[df1[col_name].str.contains(ch_kw2)]
df206 = pd.concat([df204,df205])
print('該当するデータ件数 = ',len(df206))

特定の数値項目で条件検索する
今度は数値データが入った列について、ある条件を満たす行を取り出すことを考えます。例えば値がゼロより大きい行を抽出するには、df[df['列名']>0]とする書き方と、df.query('列名>0')という書き方があります。

# 数値変数で条件検索する
ch_col='reviewCount'; s_val=500
df31 = df1[df1[ch_col] >= s_val]
dfi.export(df31,'top_review_no.png')
#数値検索の別の記述法
ch_col='reviewCount'; s_val=500
df32 = df1.query('{0} >= @s_val'.format(ch_col))
dfi.export(df32,'top_review_no2.png')

結果は以下のようになります。
f:id:nicjps230:20210719000425p:plain
特定の列項目が空でないケースを抽出する
次に、テキストからなる特定の列項目が空でないケースを抽出する方法を考えます。テキスト情報で何も書かれていないのは空文字列('')なので、まず空文字列を欠損値(NaN)に置き換えます。その後で欠損値を含む行を削除します。

# 列名を指定する
col_name = 'contents'
# データフレームをコピーする
df211 = df1.copy()
#空白をまずNaNに置き換え
df211[col_name].replace('', np.nan, inplace=True)
#Nanを削除 inplace=Trueでdfが上書きされる。
df211.dropna(subset=[col_name], inplace=True)
print('該当するデータ件数 = ',len(df211))

(下巻へ続く)

機械学習モデルの数学的記述と計算(その1):バイナリ学習

機械学習は、データからそれらが生成されている法則性を学んで、予測や分類に役立てるものである。
ここでは、スタンフォード大学のAndrew Ng教授の講義ノートをベースに、機械学習の数学的バックグラウンドを整理する。

教師ありバイナリ学習のデータの表現

教師ありバイナリ学習に用いられるデータを、特徴量とターゲットのペアとして考える。これを、(x,y) により表現する。ここで、 x \in \mathbb{R}^{n_{x}} は特徴量のベクトルを表し、 y \in \{0,1\} はターゲット変数を表す。
m個の学習用データの組は以下で表される


(x^{(1)},y^{(1)}),(x^{(2)},y^{(2)}),\cdots,(x^{(m)},y^{(m)})

特徴量とターゲット変数をそれぞれまとめて以下のように記述する


x=(x^{(1)},x^{(2)},\cdots,x^{(m)}),\ \ x \in \mathbb{R}^{n_{x} \times m}

y=(y^{(1)},y^{(2)},\cdots,y^{(m)}),\ \ y \in \mathbb{R}^{1 \times m}

ロジスティック回帰モデル
x \in \mathbb{R}^{n_{x}} が与えられたとき、 以下の式により \hat{y}=P(y=1|X),0 \le \hat{y} \le 1を得たいとする。
(パラメータとアウトプット)


w \in \mathbb{R}^{n_{x}},\ \ \ b \in \mathbb{R}

\hat{y}=\sigma(w ^{T} x+b)
ここで、

\sigma (z)=\frac{1}{1+\exp (-z)}

ロジスティック回帰モデルの損失関数
上記の定式化で、\{(x^{(1)},y^{(1)}),\cdots,(x^{(m)},y^{(m)})\} が与えられたとき、 \hat{y}^{(i)} \approx y^{(i)} を得たいとする。
損失(誤差)関数を以下のような「交差エントロピー誤差関数」として定義する。


J(w, b)=\frac{1}{m}\sum_{i=1}^{m}\mathcal{L}(\hat{y},y)=-\frac{1}{m}\sum_{i=1}^{m}(y^{(i)}\log \hat{y}^{(i)}+(1-y^{(i)})\log (1-\hat{y}^{(i)}))

最急降下法
上の式の J(w,b)を最小化する w,bを見つけたいとする。例えば J wの関数としてみた場合、 Jの勾配を計算して、勾配が急な方向に沿って進む。つまり以下を反復計算する


w:=w-\alpha \frac{dJ(w)}{dw}
ここで \alphaは学習率(learning rate)である。

計算グラフ
機械学習では以下のような計算グラフを考える。そして、その微分を逆伝搬で考える。これにより計算を局所化し、計算効率が高まる。
ここでは、 J(a,b,c)=3(a+bc) を例として考える。
f:id:nicjps230:20210711130639p:plain

上のロジスティック回帰モデルの計算グラフとその微分(概要のみ)は以下のようになる。
f:id:nicjps230:20210716140912p:plain

m個のデータ例でのロジスティック回帰モデル
損失関数の式を再掲する


J(w, b)=\frac{1}{m}\sum_{i=1}^{m}\mathcal{L}(a^{(i)},y^{(i)})
ここで、

a^{(i)}=\hat{y}^{(i)}=\sigma (w^{T}x^{(i)}+b)
いま、パラメータが w_{1},w_{2},bのみであるとする。上記の損失関数の、パラメータ w_{1}に関する傾きを求めるためには、以下を計算する。

\frac{\partial }{\partial w_{1}}J(w,b)=\frac{1}{m}\sum_{i=1}^{m}\frac{\partial }{\partial w_{1}}\mathcal{L}(a^{(i)},y^{(i)})

上記のように微分を逆伝搬で考える場合、各ノードの変数が1単位変化した時に最終出力変数がどのくらい変化するか、すなわち

\frac{\partial (final output var)}{\partial (var)}=dvar
を考えることになる。例えば、 w_{1}についての勾配は、上の図の微分の逆伝搬を伝っていき、

dw_{1}=\frac{\partial \mathcal{L}}{\partial w_{1}}=\frac{\partial \mathcal{L}}{\partial a}\cdot \frac{\partial a}{\partial z}\cdot x_{1}
と求められる。いま、

\mathcal{L}(a,y)=-(y\log a+(1-y)\log (1-a))
 aについての「勾配」を求めると、

da=\frac{\partial \mathcal{L}(a,y)}{\partial a}=-\frac{y}{a}+\frac{1-y}{1-a}
となる。
また、

a=\sigma (z)=\frac{1}{1+\exp (-z)}
 zについて解くと、

z=-\log (\frac{1}{a}-1)
となり、逆関数微分法より、

\frac{\partial a}{\partial z}=\frac{1}{\frac{\partial z}{\partial a}}=a(1-a)
となる。それゆえ、

dz=\frac{\partial \mathcal{L}}{\partial z}=\frac{\partial \mathcal{L}}{\partial a} \cdot \frac{\partial a}{\partial z}=(-\frac{y}{a}+\frac{1-y}{1-a}) \cdot a(1-a)=a-y
となる。
さらに、

z=w_{1}x_{1}+w_{2}x_{2}+b
より、

dw_{1}=\frac{\partial \mathcal{L}}{\partial w_{1}}=\frac{\partial \mathcal{L}}{\partial z} \cdot \frac{\partial z}{\partial w_{1}}=dz \cdot x_{1}
となる。同様に、 dw_{2}=x_{2} \cdot dz,   db=dzとなる。

アルゴリズムpythonによる実装
上記のアルゴリズムでロジスティック回帰の係数を求める計算をpythonで実装した。データはKaggle上にあるロジスティック回帰用のデータを用いた。ここでの留意点は、特徴量の線型結合の絶対値が大きくなる(数十程度)と、シグモイド関数極値に振れてしまい、対数が計算できなくなることである。そのため、ここではあらかじめダミー変数以外の変数を平均0、標準偏差1で基準化した。

# ライブラリのインポート
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import japanize_matplotlib # matplotlib日本語
# データの読み込み
df01 = pd.read_csv('/content/drive/My Drive/Colab Notebooks/Social_Network_Ads.csv')
# 性別に関するダミー変数の作成
df02 = pd.get_dummies(df01, drop_first=True, columns=['Gender'])
# ダミー変数以外の特徴量を標準化する
df02 = df02.assign(salary_std = (df02['EstimatedSalary']-df02['EstimatedSalary'].mean())/df02['EstimatedSalary'].std())
df02 = df02.assign(age_std = (df02['Age']-df02['Age'].mean())/df02['Age'].std())
# 分析に用いるデータをNumpy配列に書き出す
features_list = ['salary_std','Gender_Male','age_std']
X_train = df02[features_list].to_numpy()
y_train = df02['Purchased'].to_numpy()
# バイナリ学習のアルゴリズム
w1 = 1; w2 = 1; w3 = 1; b = 1 # 回帰係数の初期化
alpha = 1 # 学習率
m = len(y_train)
z = np.zeros(m); a=np.zeros(m); dz=np.zeros(m)
for j in range(250): # 指定した回数の反復計算
  J = 0; dw1 = 0; dw2 = 0; dw3 = 0; db = 0 # 損失関数と勾配の初期化
  for i in range(m): # データ行についての処理
    z[i] = w1 * X_train[i,0] + w2 * X_train[i,1] + w3 * X_train[i,2] + b
    a[i] = 1 / (1+np.exp(-z[i]))
    J += -(y_train[i] * np.log(a[i]) + (1-y_train[i]) * np.log(1-a[i]))
    dz[i] = a[i] - y_train[i]
    dw1 += X_train[i,0] * dz[i]
    dw2 += X_train[i,1] * dz[i]
    dw3 += X_train[i,2] * dz[i]
    db += dz[i]
  J /= m; dw1 /= m; dw2 /= m; dw3 /= m; db /= m
  w1 -= alpha * dw1
  w2 -= alpha * dw2
  w3 -= alpha * dw3
  b -= alpha * db
  if j % 20 == 0:
    print('j={}, J={:.3g}, w1={:.3g}, w2={:.3g}, w3={:.3g},b={:.3g}'.format(j,J,w1,w2,w3,b))

損失関数Jの収束の様子から学習率αを調節していく。αを調節した結果は以下のようになった。収束には100回以上の反復計算を要している。
f:id:nicjps230:20210715140252p:plain