介護離職を予測してみた



# 目次
[1. はじめに](#はじめに)
[2. 本記事の概要](#本記事の概要)
[3. 分析開始まで](#分析開始まで)
[4. 作成したプログラム](#作成したプログラム)
[5. 最終的な精度評価](#最終的な精度評価)
[6. 考察](#考察)
[7. 今後の活用・展望](#今後の活用・展望)


## はじめに

ご覧くださりありがとうございます。
Aidemy Premiumにて「データ分析コース」を受講いたしました。
2年ほど前にG検定を取得したことから、本業ではあまり使わないものの、統計やAI分野に興味を持っており、まとまった時間が出来たことから受講をしてみることにしました。

最終成果物の制作が修了要件とのことで、やるのであれば在り来りなものよりは少しは独自性があるものにチャレンジしたいと思いました。そこで、副業を通して持っていた介護に関する調査データを活用することにしました。
※データ利用については副業先から許可を得ております。

**簡単な自己紹介**
大学で統計学をかじったアラサー
本業は事業企画なのでデータは関係あったりなかったり程度
Excelでも十分担当業務は完結できる

＜　所有資格　＞
- 統計検定2級
- 基本情報技術者
- G検定

**なお、「このブログはAidemy Premiumのカリキュラムの一環で、受講修了条件を満たすために公開しています」**


## 本記事の概要

初めに筆者は統計学を少しかじった程度であり、AIは素人ですので粗い点の多いことにはご了承頂ければ幸いです。
- 対象読者：
    - これからAidemyなどを通してAIの学習をしてみたい方
    - 介護に近しい分野で関心を持って頂ける方
    ➔　介護関連のデータ活用はこれからも私のテーマなので
- この記事で分かること
    - Aidemyデータ分析コースの感想・レビュー
    - 初学者が真面目に課題と向き合ったときの出来栄え（限界感じました笑）
- この記事で扱わないこと
    - オープンソースデータの活用

## 分析開始まで
### コースの進捗
9/24 ：受講開始
10/8 ：データ分析コース成果物課題を残して学科受講が完了
10/15：成果物課題　初回カウンセリング
11/08：成果物課題　提出

### 個人的に伝えたいポイント
データ分析コースの学科履修はサクサク進めることができ、最終成果物の制作までに用した時間は40時間弱でした。シラバスに記載されている時間の半分くらいで進んできたので、少々拍子抜けしながら駆け抜けました。
成果物に取り組みながら思ったのですが、成果物を作りながら講座を受けたほうが実りがあったように感じました。もちろん、最初はコースの概要も掴みきれない中で進むので無理もないのですが、結局成果物を作りながらGoogle検索、ChatGPTを活用しました。Aidemyのコースを見返すことも多々ありました。

### 分析方針
<b>最終成果物については、実践的な内容でやってみたい気持ちがあり、数年前に取得したきりなんとなくもったいないなと思っていた、「介護離職」調査のデータを活用してみたいと思い当たりました。</b>

離職有無を表すカラムがあったので、これをy_labelとして使い、**回帰（分類）** のモデルを作ることにしました。
X_labelについては活用するデータを取得した調査を実施した際には回答者の介護に対する考え方なども聴取していましたが、実際にそのような指標を現実世界で聴取できるかと言えばそうではないので使わないことにしました。
ただし、データを削ったとしてもパラメータが51残っており、先が思いやられましたが、とりあえずやってみることにしました。

### 開発環境
Surface PRO 7
Microsoft Windows 11 Home
Google Colaboratory
Python 3.10.12

### 使用データ
**サンプル数**
1600（内訳：各400）
　　├　①現在介護をしている介護者かつ、介護離職あり
　　├　②現在介護をしている介護者かつ、介護離職なし
　　├　③過去に介護をしていた介護者かつ、介護離職あり
　　└　④過去に介護をしていた介護者かつ、介護離職なし

**カラム**
<details><summary>クリックで表示します</summary>
データカラム
SAMPLEID<br>
SEX<br>
AGE<br>
AGEID<br>
PREFECTURE<br>
AREA<br>
MARRIED<br>
CHILD<br>
HINCOME<br>
PINCOME<br>
CELL<br>
CELLNAME<br>
PRESENT JOB<br>
CARE RECEIVER<br>
NUM EMPLOYEE<br>
INDUSTRY<br>
OCCUPATION<br>
LABOR<br>
WORKING TIME<br>
介護休業(長期)<br>
介護休暇(短期)<br>
SHORT TIME<br>
FLEX<br>
REMOTE<br>
CARE LEVEL START<br>
CARE LEVEL<br>
DEMENTIA<br>
CR AGE<br>
CR SEX<br>
CARE EXPENSE<br>
CARE SITUATION<br>
TIME REQUIRED<br>
CARE FREQ<br>
TIME SPENT<br>
地域包括 活用<br>
COHABITANT<br>
COHAB YOUNG<br>
BUSINESS SCALE<br>
CARE YEAR<br>
同居別居<br>
制度利用有無<br>
介護開始<br>
介護関わり始め<br>
介護関わり終わり<br>
介護負担始め<br>
介護負担終わり<br>
介護離職<br>
現在開始<br>
過去初め<br>
過去終わり<br>
離職時期<br>
話し相手<br>
買い物<br>
洗濯<br>
掃除<br>
散歩<br>
服薬の手助け<br>
食事介助<br>
調理を含む食事の準備・後始末<br>
排せつ介助<br>
体位交換・起居<br>
入浴介助<br>
着替え<br>
洗髪<br>
CARE RELATION_SPOUSE<br>
CARE RELATION_PARENT<br>
CARE RELATION_SPARENT<br>
CARE RELATION_GRANDPARENT<br>
CARE RELATION_SGRANDPARENT<br>
CARE RELATION_CHILD<br>
</details>

## 作成したプログラム
### データの取得
Google Driveからデータを取得していきます。
```python:
from google.colab import drive
drive.mount('/content/drive')
df = pd.read_csv("/content/drive/MyDrive/004_Data Science engineer/Aidemy/data analysis final output/data_aidemy3.csv")
```

### 活用するライブラリ一覧
ここに記載したのは最終的に活用することになったライブラリです。だんだん増えていったので最終的にまとめました。
```python:
import matplotlib.pyplot as plt
import seaborn as sns
import numpy as np
import pandas as pd
import xgboost as xgb
from sklearn.model_selection import train_test_split
from sklearn.metrics import f1_score, confusion_matrix
from sklearn.linear_model import LogisticRegression, LinearRegression
from sklearn.svm import LinearSVC, SVC
from sklearn.tree import DecisionTreeClassifier
from sklearn.ensemble import RandomForestClassifier
from sklearn.neighbors import KNeighborsClassifier
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LinearRegression
from sklearn.svm import SVR
from sklearn.tree import DecisionTreeRegressor
from sklearn.ensemble import RandomForestRegressor
from sklearn.metrics import mean_squared_error
from sklearn.metrics import log_loss
from sklearn.metrics import mean_squared_error

```

### データの中身確認
~~正直Excelのほうが早くないかと思いつつ、我慢してPythonでやってみました。~~
今回はカラム数が多く、デフォルト設定だとdescribeの出力結果が短縮表記されてしまいます。それを解消するために最初の1行をはさみました。
```python:
pd.set_option('display.max_columns', 100)
df.describe()
df.head()
```
結果、知ってはいましたが現在介護者と過去介護者で聴取項目が違っていたり（＝カラムによっては半分が欠損値になっている）、いらないカラムがあったり、漢字で入力されたカラムがあったり、様々データのクレンジング、ハンドリングが必要になることがわかりました。

パラメータの選別をしたほうがいいかと思っていたのですが、AidemyのSlack質問を使って聞いてみたところ、
>まずは全部使ってみて、精度がうまく出なかった場合に特徴量を減らしていってみましょう。特徴量の減らし方ですが、基本的にはモデルがその目的変数を予測するのに重要な項目をピックアップして残すのが一般的です。
ピックアップの仕方は普段介護に関連したお仕事やご活動をされているなら普段活動している中でこの項目を大事だと気づかれていることがありましたらそういったものを残したり、もっと機械的な方法としては目的変数との相関を調べて相関が高い特徴量を残すなどの方法もあります。

とのお言葉をいただき、全て突っ込んでみることにしました。

### 変数のチェックと変換（①）
チューターさんの言葉通り捉え欠損値は0で埋めようという非常に安易な考えでいました。
そこで、最低限の変数の変換を行うことにします。
```python:
#離職したかどうかのカラムを作りたいので、上記カテゴリ①～④で表したカラム【CELL】を変換する
df["CELL"].replace({2: 0, 4: 0, 3: 1}, inplace=True)
#カラム【同居別居】について、バイナリに変換
df["同居別居"].replace({"同居": 1, "別居": 0}, inplace=True)
#カラム【MARRIED】について、バイナリに変換
df["MARRIED"].replace({1: 0, 2: 1}, inplace=True)
```
### モデリング（①）
明らかに必要なさそうなカラム、被っているカラムを消してモデリングを行うデータフレームとする
```python
df_1 = df.drop(['PRESENT JOB','PREFECTURE','AREA','AGEID','CELLNAME', 'SAMPLEID',  '介護開始', '介護関わり始め', '介護関わり終わり', '介護負担始め', '介護負担終わり', '現在開始', '離職時期', '介護離職', '過去初め', '過去終わり','CARE LEVEL START','BUSINESS SCALE'], axis =1)
```

モデルはAidemy内の講座で出てきた5つの手段を試しました。
```python
#　とりあえず全カラムぶちこみモデル化
#df_1は欠損値を0で補完する
df_1.fillna(0, inplace=True)
#これにてCELLをy_labelとしてセットする準備が整った
y = df_1.CELL
X = df_1.drop(columns=['CELL'])
# データを学習に使う分と評価の分に分ける
train_X, test_X, train_y, test_y = train_test_split(X, y, test_size=0.3, random_state=42)

def train_and_evaluate(model, train_X, train_y, test_X, test_y):
    model.fit(train_X, train_y)
    pred_y = model.predict(test_X)
    f1 = f1_score(test_y, pred_y)
    conf_matrix = confusion_matrix(test_y, pred_y)

    print(f"F値: {f1}")
    print("混同行列:\n", conf_matrix)

# ロジスティック回帰
model_lr = LogisticRegression(max_iter=1000)
train_and_evaluate(model_lr, train_X, train_y, test_X, test_y)

# 線形SVM
model_lsvm = LinearSVC(max_iter=1000)
train_and_evaluate(model_lsvm, train_X, train_y, test_X, test_y)

# 非線形SVM
model_svm = SVC()
train_and_evaluate(model_svm, train_X, train_y, test_X, test_y)

# 決定木
model_dtc = DecisionTreeClassifier()
train_and_evaluate(model_dtc, train_X, train_y, test_X, test_y)

# ランダムフォレスト
model_rfc = RandomForestClassifier(max_depth=11,n_estimators=38,random_state =14)
train_and_evaluate(model_rfc, train_X, train_y, test_X, test_y)

# k-NN
model_knn = KNeighborsClassifier()
train_and_evaluate(model_knn, train_X, train_y, test_X, test_y)
```
**出力された結果**
この時点でランダムフォレストで相対的にに良い結果が出た。
ここからは欠損値の保管や、カテゴリ変数の変換などの所作を行い、モデル精度を高めようと考えた。
```
ロジスティック回帰のF値: 0.6487603305785123
ロジスティック回帰の混同行列:
[[153  85]
 [ 85 157]]
線形SVMのF値: 0.6142857142857143
線形SVMの混同行列:
[[189  49]
 [113 129]]
非線形SVMのF値: 0.0
非線形SVMの混同行列:
[[238   0]
 [242   0]]
決定木のF値: 0.5833333333333334
決定木の混同行列:
[[140  98]
 [102 140]]
ランダムフォレストのF値: 0.6983471074380165
ランダムフォレストの混同行列:
[[165  73]
 [ 73 169]]
k-NNのF値: 0.5828092243186583
k-NNの混同行列:
[[142  96]
 [103 139]]
 ```

### 欠損値の対応
```python:
#欠損値等チェック
pd.set_option('display.max_columns', 100)
df.describe()
```
以下のカラムで欠損値が確認できた。
|カラム|情報|状況把握|取った対応|
|---|---|---|---|
|DEMENTIA|認知症有無|不明など未回答者あり|(A)予測により補完|
|HINCOME|世帯年収|未回答者あり|(A)予測により補完|
|PINCOME|個人年収|未回答者あり|(A)予測により補完|
|TIME REQUIRED|要介護者までの所要時間|同居介護者は未回答|0で補完|
|VISIT FREQ|介護訪問回数|通い介護者のみ回答|(B)CARE FREQと同義なので合体する|
|CARE FREQ|介護回数|同居介護者のみ回答|(B)VISIT FREQと同義なので合体する|
|TIME SPENT|1日平均介護時間|介護に携わっていない場合は未回答|0で補完|

以下、(A)(B)のコードについて記載
### (A)への対応HINCOMEの例
```python:
df_dem = df.copy()
df_dem = df_dem.drop(['PRESENT JOB', 'PREFECTURE', 'AREA', 'AGEID', 'CELLNAME', 'SAMPLEID', '介護開始', '介護関わり始め', '介護関わり終わり', '介護負担始め', '介護負担終わり', '現在開始', '離職時期', '介護離職', '過去初め', '過去終わり', '制度利用有無', 'CARE LEVEL START', 'BUSINESS SCALE'], axis=1)

# 準備
X_train_dem = df_dem.dropna(subset=['HINCOME'])  # 欠損値のない行
X_test_dem = df_dem[df_dem['HINCOME'].isnull()]  # カラムHINCOMEが欠損値の行

# とりあえず欠損値に0を入れる
X_train_dem.fillna(0, inplace=True)
X_test_dem.fillna(0, inplace=True)

# モデルの選択 (ここではランダムフォレストを使用)
model_dem = DecisionTreeClassifier()

# モデルのトレーニング
model_dem.fit(X_train_dem.drop(['HINCOME'], axis=1), X_train_dem['HINCOME'])

# 予測
predicted_dem = model_dem.predict(X_test_dem.drop(['HINCOME'], axis=1))

# インデックスをリセットする
df_dem.reset_index(drop=True, inplace=True)
# 'HINCOME' カラムの値が 0 である行のインデックスを取得
zero_income_indices = df_dem[df_dem['HINCOME'] == 0].index
# 予測値を代入
for i, idx in enumerate(zero_income_indices):
  df_dem.at[idx, 'HINCOME'] = predicted_dem[i]

# データフレーム `df` に結果を反映
df['HINCOME'] = df_dem['HINCOME']
value_counts_cell = df['HINCOME'].value_counts()
value_counts_cell_df = pd.DataFrame({'Value': value_counts_cell.index, 'Count': value_counts_cell.values})
# データフレームを表示
print(value_counts_cell_df)
```

### カテゴリカル変数の対応
one-hot encodingを用いて対応しました。以下はその一例です。
```python:
fig_sex = sns.countplot(x='SEX', hue='CELL', data = df)
plt.show()
# SEXカラムをOne-Hotエンコーディング
df_2 = pd.get_dummies(df_2, columns=['SEX'], prefix=['SEX'],drop_first=True)
```
one-hot encodingでは、生成するカラム数をカテゴリ-1にする必要があるということも調べる中でわかりました。その所作も上記コードの中で行っています。

### グラフ化 → カラムの作成
各カラムをグラフ化する中で傾向が見られたものについては、その傾向を特徴量として表現してみました。
以下はその一例としてHINCOME（世帯年収）に関するコードとグラフを示します。
```python:
fig_pref = sns.countplot(x='HINCOME',hue='CELL', data = df)
plt.show()
```
**出力されたグラフ**
![P,HINCOME可視化.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3577454/db6e946e-ba38-e071-5f2e-7640971fd8e9.png)

上記グラフからは1～3のときに、離職者が比較的多く存在し、4～7のときに未離職者が比較的多く存在していそうだと見て取ることができました。
そこで以下のように特徴量を増やしました。
```python:
#HINCOME
df_2['Low_HINCOME'] = df_2['HINCOME'].map(lambda s: 1 if 1 <= s <= 3 else 0)
df_2['Med_HINCOME'] = df_2['HINCOME'].map(lambda s: 1 if 4 <= s <= 7 else 0)
```

### 標準化
モデル制作を進めていると、SVMなど距離をベースとするモデルに関しては、スケーリングをするほうがベターであるということがわかりました。
以下はAGE＝年齢について、標準化をしています。
```python:
scaler = StandardScaler()
df_2['AGE'] = scaler.fit_transform(df_2[['AGE']])
```


## 最終的な精度評価
以上で紹介してきた所作を各カラムに対して適用し、再度モデリングを実施しました。ほんの僅かではありますが、最初に実施したモデリングのF値と比較して精度が良くなったようです。
```
ロジスティック回帰のF値: 0.6467661691542289
ロジスティック回帰の混同行列:
 [[128  71]
 [ 71 130]]
線形SVMのF値: 0.6773455377574371
線形SVMの混同行列:
 [[111  88]
 [ 53 148]]
非線形SVMのF値: 0.0
非線形SVMの混同行列:
 [[199   0]
 [201   0]]
決定木のF値: 0.6122448979591837
決定木の混同行列:
 [[128  71]
 [ 81 120]]
ランダムフォレストのF値: 0.7067669172932332
ランダムフォレストの混同行列:
 [[142  57]
 [ 60 141]]
k-NNのF値: 0.6005089058524175
k-NNの混同行列:
 [[125  74]
 [ 83 118]]
```

## 考察
最後に考察を述べていきたいと思います。
まず介護離職の予測を、一般的に行政窓口でも把握できるようなデモグラフィックな情報から試みました。結果として一定の予測精度を実現することができました。
介護に潜むリスクは一時的なデータから見積もれるほど簡単ではないですが、データによってその裏付けや、検討の補助ができることは意味があると考えます。

一方で、今回は様々にデータハンドリングを試みることで予測精度の向上を目指しましたが、ほとんど精度を上げることができず悔いが残りました。偶に耳にする、データ分析におけるデータ収集、データの質が重要であるということを身をもって経験することになったと感じています。
できればF値を0.8くらいにまで引き上げたく、様々な試行錯誤をしましたが、同頑張っても0.71を上回ることはありませんでした。

しかし、今回モデリングを反復する中でRandom Forestにおける各説明変数の重要度を調べる方法を知ることもできました。参考までにそのコードを記載します。
```python:
fi = model_rfc_2.feature_importances_
idx = np.argsort(fi)[::-1]
top_cols, top_importances = train_X2.columns.values[idx][:50], fi[idx][:50]
# カラムと重要度をデータフレームにまとめる
feature_importance_df = pd.DataFrame({'Feature': top_cols, 'Importance': top_importances})
# すべてのカラムを表示
pd.set_option('display.max_rows', None)
print('Random Forest Importance for All Columns')
print(feature_importance_df)
```
上記のコードから、今回の予測モデル構築に当たって、どの説明変数が重要だったかを可視化することができました。介護離職に対する重要度は経験があったとしても複合的で判断が難しいため、こういった分析も役に立てることができそうだと感じました。モデリングを反復する中で上記のような発見ができることも、データ分析の奥深さだと感じました。


## 今後の活用・展望
Aidemyでの学習は実はこれからも続きます。私は６ヶ月コースを受講することで、学び放題プランというオトクなプランを使っており、E検定に合格するところまでをゴールにしているためです。なので短期的な目標はE検定を2月に受け、合格することが目標となります。
もう少し先に向けては、やはりデータ分析は面白いですし、持ち腐れにならないように案件をもらえるような能力に磨いていきたいと思います。