メモ(中島)
の編集
Top
/
メモ(中島)
[
トップ
] [
編集
|
差分
|
履歴
|
添付
|
リロード
] [
新規
|
一覧
|
検索
|
最終更新
|
ヘルプ
]
-- 雛形とするページ --
島部/menu/template
[[中島]] ---- ~ &size(15){''[[参考文献とかメモとか>./参考文献とかメモとか]]''};~ ~ ---- #memo([[技術資料]]\n*目次 [#f3c34276]\n\n#CONTENTS\n\n*本研究の目的 [#cec342b3]\n本研究の目的は，ヘドニック・アプローチにおける限界を克服し，価格形成におけるより\n正確で信頼性の高い評価手法を提供することである．ヘドニックモデルでは，多重共線性，\n欠落変数のバイアス，交互作用効果および非線形効果の取り扱いに問題がある[2]．これら\nは不動産市場や消費財の価格評価において誤差を生じさせる要因となる．本研究は，これ\nらの問題を解決するために，スパース推定法を活用した新たなアプローチを提案する．\n従来のヘドニックモデルは，そのシンプルな線形構造と主要変数の関係を明確に示す点\nで強力であったが，現実の市場では非線形な相互作用や高次の効果が絡むため，これらを\n十分に捉えることができないことが多かった．また，異なる要因が相互作用を持つ場合，そ\nの効果を正確に評価するためには，より複雑なモデルが求められる．そこで本研究では，\nElasticNetおよびその進化形であるAdaptiveElasticNet(AEN)を用い，従来の限界を克服\nすることを目指す．\nElasticNetはL1およびL2正則化を組み合わせることで，多重共線性の問題に対応する\nと同時に，重要な変数の選択を行う能力を持つ．これにより，データ中で重要な特徴を捉\nえ，効果的な変数選択を実現できる．また，AENは，従来のElasticNetをさらに発展させ，\nデータの特性に応じて正則化パラメータを動的に調整する．この適応的なアプローチは，複\n雑なデータ構造に対して非常に効果的であり，非線形性や交互作用を持つデータにも適用\nできる．これにより，不動産市場や消費財の価格形成における微細な変動や相互作用を捉\nえ，価格予測の精度を向上させることが可能となる．\nさらに，本研究では，従来の線形回帰モデルでは捉えきれなかった非線形性や交互作用\n効果をスパース推定法を活用することで反映させることができる．非線形回帰や機械学習\n技術と組み合わせることにより，より精緻な価格形成モデルが構築され，従来の方法では\n見過ごされがちな要因を捉えることが可能となる．例えば，不動産市場における立地や面\n積の相互作用，消費財におけるブランド力と機能性の複合的影響など，これらの複雑な関\n係性を反映させることで，より現実的な市場分析が行えるようになる．\nまた，こうした高度な統計手法や機械学習技術の活用により，従来のヘドニック・アプ\nローチを大きく強化することができる．特に，従来のアプローチでは捉えきれなかった価\n格形成における微細な要因を，スパース推定を用いることで明示化し，より高精度な価格\n評価を可能にする．これにより，不動産市場や消費財市場における価格形成メカニズムの\n理解が深まり，実務的な応用に貢献することが期待される．\n*使用するファイル [#pe64d1b6]\n本研究で使用するファイルはすべてGoogleDriveにある。\n&br;\nuser:iie.lab.tpu.2324@gmail.com\n&br;\n\n「マイドライブ」→「学生」→「24_o3中島」→「保存期間5年」→「３）卒論（プログラム）」の中にある。\n&br;\n「引継ぎ用.zip」をダウンロードする。\n&br;\nファイルの中身は以下のようになっている。\nこれを上から順に実行していく過程を以下の章に示す。\n\n|フォルダ名|用途|ファイル名|ファイルの場所|\n|データセットの作成|犯罪データの取得|get_crime_data.py|データセットの作成\get_crime_data.py|\n|データセットの作成|統計データなどの作成|hop_統計.ipynb|データセットの作成\hop_統計.ipynb|\n|データセットの作成|取引データのメッシュ変換|緯度経度.py|データセットの作成\緯度経度.py|\n|データセットの作成|各データを取引データにまとめる|hop_メッシュ.ipynb|データセットの作成\hop_メッシュ.ipynb|\n|データセットの作成|犯罪場所との最短距離|distance.ipynb|データセットの作成\distance.ipynb|\n|前処理|データセットの前処理|前処理.ipynb|前処理\前処理.ipynb|\n|変数選択|変数の選択とハイパーパラメータのチューニング|変数選択__削除.ipynb|変数選択\変数選択_削除.ipynb|\n|要因分析|選択された変数を用いたヘドニック・アプローチによる要因分析|hedo_削除.ipynb|要因分析\hedo_削除.ipynb|\n~\nコードの変更点はディレクトリを変更するだけで大丈夫.\n\npythonのバージョンは12.0.0を用いる.\n\n使用したモジュールのバージョンは下になります。（すべては使用しておりません）\n#ref(モジュール１.png,,)\n#ref(モジュール2.png,,)\n#ref(モジュール3.png,,)\n#ref(モジュール4.png,,)\n** 犯罪発生データの取得 [#k37ae06a]\n\n「''get_crime_data.py''」を実行すると、富山県警が公開している「犯罪発生マップ」（http://www.machi-info.jp/machikado/police_pref_toyama/index.jsp）から犯罪発生データを取得し、「''crimes.csv''」として保存されます。\n&br;\n\n**建築年数の算出 [#i8ed6a27]\n建築年の年号を西暦にするために以下の式を書いた列を挿入する。\n&br;\n =IF(LEFT(M2,2)="令和", 2018 + MID(M2,3,LEN(M2)-3),\n IF(LEFT(M2,2)="平成", 1988 + MID(M2,3,LEN(M2)-3),\n IF(LEFT(M2,2)="昭和", 1925 + MID(M2,3,LEN(M2)-3), "エラー")))\n次に取引時点から西暦を抜き出す。\n&br;\n =LEFT(W2,4)\nそして、取引時点の西暦から建築年の西暦を引くことで建築年数を算出する。\n\n**使用するデータ\n実際の取引データはhttps://www.reinfolib.mlit.go.jp/realEstatePrices/ このサイトからダウンロードする。最新版を用いる場合はこちらからダウンロードを。\n\nデータセットの作成\data\estats_mesh、データセットの作成\data\estats_blockのデータはhttps://www.e-stat.go.jp/から取得できる。最新版を取得する場合はこちらから。\nhttps://www.e-stat.go.jp/gis/statmap-search?\npage=1&type=1&toukeiCode=00200521&toukeiYear=2020&aggregateUnit=H&serveyId=H002005112020&statsId=T001108&datum=2000&prefCode=16\nhttps://www.e-stat.go.jp/gis/statmap-search?page=1&type=1&toukeiCode=00200521&toukeiYear=2020&aggregateUnit=H&serveyId=H002005112020&statsId=T001101&datum=2000&prefCode=16 この二つのurlからメッシュごとのデータ。\nhttps://www.e-stat.go.jp/stat-search/files?page=1&layout=datalist&toukei=00200521&tstat=000001136464&cycle=0&tclass1=000001136472&tclass2=000001159889&cycle_facet=tclass1%3Acycle&tclass3val=0 このurlから章地域ごとのデータを。\n** データの作成と収集（統計データや施設データ）[#x2b1cdcf]\n「''hop_統計.ipynb''」を実行すると、取引データのメッシュごとの統計データの作成が行われる.\nNAVITIMEからのスクレイピングや地図画像からのデータ作成を行い、メッシュごとのデータを最終的にdata_by_mesh.csvに格納する.\n\n&br;\n df_torihiki = pd.read_csv(\n HERE / "df_housing_toyamacityメッシュ入り.csv",\n index_col = "住所",\n encoding = 'cp932') \n\nメッシュコードを文字列として格納 \n&br;\n meshes = df_torihiki["mesh_code"].astype(str).tolist()\n作成する統計データは以下のようになる。\n#ref(作成する統計データ.png,,)\n&br;\n import pandas as pd\n # ここで上で指定した説明変数を生成\n\n df_estats_block = pd.read_csv(\n DATA_DIR / "estats_block.csv",\n index_col = "KEY_CODE")\n df_estats_mesh = pd.read_csv(\n DATA_DIR / "estats_mesh.csv",\n index_col = "KEY_CODE")\n\n df_selected_block = pd.DataFrame()\n df_selected_mesh = pd.DataFrame()\n\n for line in stats_data_text.splitlines():\n \n formula, label = line.rsplit(", ", 2)\n \n if len( formula.split(" ") ) == 1:\n \n if formula in df_estats_block.columns:\n \n df_selected_block.loc[:, label] = df_estats_block.loc[:, formula]\n \n if formula in df_estats_mesh.columns:\n \n df_selected_mesh.loc[:, label] = df_estats_mesh.loc[:, formula]\n \n elif len( formula.split(" ") ) >= 2:\n \n contains_keys = re.findall("T.........", formula)\n \n if contains_keys[0] in df_estats_block.columns:\n \n for key in contains_keys:\n \n formula = formula[:formula.find(key)] + "df_estats_block.loc[:, '" + \\n formula[formula.find(key) : formula.find(key) + len(key)] + \\n "']" + formula[formula.find(key) + len(key):]\n \n exec(f"df_selected_block.loc[:, '{label}'] = {formula}")\n \n elif contains_keys[0] in df_estats_mesh.columns:\n \n for key in contains_keys:\n \n formula = formula[:formula.find(key)] + "df_estats_mesh.loc[:, '" + \\n formula[formula.find(key) : formula.find(key) + len(key)] + \\n "']" + formula[formula.find(key) + len(key):]\n \n exec(f"df_selected_mesh.loc[:, '{label}'] = {formula}")\n \n df_selected_mesh = df_selected_mesh.fillna(0)\n df_selected_block = df_selected_block.fillna(0)\n df_selected_mesh\n df_selected_block\n&br;\n # 次は施設データ．\n # 統計データと同じで，カンマの左側の数字はNAVITIMEのカテゴリID\n # コンビニ一覧のURLは "https://www.navitime.co.jp/category/0201/" なので，カテゴ \n リIDは「0201」\n\n building_data_text = '''\\n 0201, コンビニ\n 0802001, 駅\n 0805, 駐車場\n 0812, 駐輪場\n 0501, 金融機関\n 0603, 旅館\n 0608, ホテル\n 0202, スーパー\n 0205, ショッピングモール\n 0204, デパート\n 0502002, 警察署／交番\n 0504001, 小学校\n 0504002, 中学校\n 0504003, 高校\n 0504004, 大学\n 0504006, 保育園／幼稚園\n 0101, レジャー施設'''\n&br;\n # ここで，NAVITIMEからスクレイピング\n\n df_buildings = pd.DataFrame()\n\n def address2coords(address: str):\n \n for trials in range(5):\n \n res = requests.get(\n urllib.parse.quote(\n string = f"https://msearch.gsi.go.jp/address-search/AddressSearch?q={address}",\n safe = "://?="\n )\n )\n \n if res.status_code == 200:\n \n break\n \n else:\n \n time.sleep(5 * trials)\n \n else:\n \n raise Exception("国土地理院のジオコーディングAPIに接続できません。")\n \n result = res.json()\n \n time.sleep(0.5)\n # サーバに負担をかけないように0.5秒待つ\n \n if len(result) == 0:\n \n return None\n \n else:\n \n return result[0]["geometry"]["coordinates"]\n # 複数ある場合でも，0番目を採用\n\n def get_building_info(category):\n \n df_buildings = pd.DataFrame(columns = ["name", "address", "lng", "lat"])\n\n for i in range(1, sys.maxsize):\n # "while True"と同義\n \n if i == 1:\n \n url = f"https://www.navitime.co.jp/category/{category}/16201/"\n \n else:\n \n url = f"https://www.navitime.co.jp/category/{category}/16201/?page={i}"\n \n for trials in range(5):\n \n res = requests.get(url)\n \n if res.status_code == 200:\n \n break\n \n else:\n \n time.sleep(5)\n \n else:\n \n break\n # NAVITIME から 5 回連続で取得できなかったら、終了\n \n soup = BeautifulSoup(res.text, "html.parser")\n \n buildings = []\n for element in soup.find_all(class_ = "spot-name-text"):\n \n content = element.get_text()\n buildings += [{"name": content}]\n \n count = 0\n for element in soup.find_all(class_ = "spot-detail-value-text"):\n \n content = element.get_text()\n \n matches = re.findall(r"(北海道|東京都|大阪府|京都府|神奈川県|和歌山県|鹿児島県|..県)(.{1,6}[市郡区町村])", content)\n # 完璧ではないが、これで住所かを判定\n \n if len(matches) > 0:\n \n address = content\n coords = address2coords(address)\n print(address, coords)\n \n buildings[count]["address"] = address\n \n if coords == None:\n buildings[count]["lng"] = None\n buildings[count]["lat"] = None\n else:\n buildings[count]["lng"] = coords[0]\n buildings[count]["lat"] = coords[1]\n \n count += 1\n \n for building in buildings:\n \n df_buildings = add_dict_to_df(building, df_buildings)\n \n return df_buildings\n\n for line in building_data_text.splitlines():\n \n code, label = line.rsplit(", ", 2)\n \n if not (DATA_DIR / f"{code}.csv").exists():\n \n df_this_buildings = get_building_info(code)\n df_this_buildings.to_csv(DATA_DIR / f"{code}.csv", index = False)\n \n else:\n \n df_this_buildings = pd.read_csv(DATA_DIR / f"{code}.csv")\n \n df_buildings[label] = None\n \n for index, row in df_this_buildings.iterrows():\n \n if not index in df_buildings.index:\n \n df_buildings.loc[index] = pd.Series(dtype='float64')\n \n df_buildings.at[index, label] = f"{row['lng']},{row['lat']}"\n\n&br;\n df_buildings\n\n&br;\n # ここでグリッドセルごとのデータを作成する\n # さっき作った統計データ：そのまま入れる\n # さっき作った施設データ：そのグリッドセルにある数と，そのグリッドセルまでの \n 最短距離を求め，入れる\n # 地図画像にもとづく特徴量：\n # MapboxのAPIを使って，地図画像を取得する．\n # 特定のRGB値（たとえば，建物だったら[255, 0, 0]）のピクセル数を求め，画像 \n サイズの512 * 512で割って（面積比率），入れる．\n # さらに，道路だけを抜き出した2値画像をつくり，道路ネットワークを抽出．ノ \n ード数，エッジ数，密度，平均次数を求めて，入れる．\n\n df_data_by_mesh = pd.DataFrame(columns = ["mesh_code"])\n df_data_by_mesh = df_data_by_mesh.set_index("mesh_code")\n\n _count = 1\n for meshcode, coords in meshes.items():\n \n outputs = get_data_in_this_mesh(meshcode, coords)\n df_data_by_mesh = add_dict_to_df(outputs["data"], df_data_by_mesh, \n outputs["meshcode"])\n \n print(f"{_count} / {len(meshes)}")\n _count += 1\n取引データに存在するメッシュごとの統計データを作成する。\n&br;\n df_data_by_mesh\n df_data_by_mesh.to_csv(HERE / "data_by_mesh.csv")\n \n**緯度経度からメッシュを格納 [#kea42b7d]\n\n「''緯度経度.py''」を実行すると、地理情報（緯度・経度）から日本独自のメッシュコードを計算し、住宅データなどのCSVファイルにそのコードを追加し、「''df_housing_toyamacityメッシュ入り.csv''」を生成します。\n&br;\n import pandas as pd\n def latlon_to_mesh(lat, lon):\n """\n 緯度経度から9桁のメッシュコード（2分の1地域メッシュ）を取得する。\n :param lat: 緯度（度）\n :param lon: 経度（度）\n :return: 9桁のメッシュコード（文字列）\n """\n # 1次メッシュ\n mesh1_lat = int(lat * 1.5)\n mesh1_lon = int(lon - 100)\n # 2次メッシュ\n mesh2_lat = int((lat * 1.5 - mesh1_lat) * 8)\n mesh2_lon = int((lon - 100 - mesh1_lon) * 8)\n # 3次メッシュ\n mesh3_lat = int(((lat * 1.5 - mesh1_lat) * 8 - mesh2_lat) * 10)\n mesh3_lon = int(((lon - 100 - mesh1_lon) * 8 - mesh2_lon) * 10)\n # 4次メッシュ（2分の1地域メッシュ）\n mesh4_lat = int((((lat * 1.5 - mesh1_lat) * 8 - mesh2_lat) * 10 - mesh3_lat) * 2)\n mesh4_lon = int((((lon - 100 - mesh1_lon) * 8 - mesh2_lon) * 10 - mesh3_lon) * 2)\n # メッシュコードの組み立て（9桁になるよう修正）\n mesh_code = f"{mesh1_lat:01}{mesh1_lon:01}{mesh2_lat:01}{mesh2_lon:01}{mesh3_lat:01}{mesh3_lon:01}{mesh4_lat:01}{mesh4_lon:01}"\n # 最後の2桁のみ取り出して、9桁にする\n mesh_code = mesh_code[:9]\n return mesh_code\n # CSVファイルの読み込み\n file_path = r"C:\Users\nt011\Desktop\研究\富山取引のほう富山市のみ \n \df_housing_toyamacity.csv"\n df = pd.read_csv(file_path, encoding='utf-8')\n # メッシュコードの計算\n df["mesh_code"] = df.apply(lambda row: latlon_to_mesh(row["緯度"], row["経度"]), axis=1)\n cols = ["mesh_code"] + [col for col in df.columns if col != "mesh_code"]\n df = df[cols]\n # 結果のCSVを保存\n output_path = r"C:\Users\nt011\Desktop\研究\富山取引のほう富山市のみ \n \df_housing_toyamacityメッシュ入り.csv"\n df.to_csv(output_path, index=False)\n print(f"処理完了: {output_path}")\n#ref(作成されるデータ.png,,)\n\n\n\n「''hop_メッシュ.ipynb''」を実行すると、メッシュコードごとにデータを結合し「''取引_統計データ付き.csv''」に生成します。緯度・経度をもとに、500mメッシュ単位で犯罪データを分類できます。\n # メッシュGeoJsonを，分かりやすい形に直す\n _meshes = {}\n for mesh in meshes:\n meshcode = mesh["properties"]["code"]\n coords = mesh["geometry"]["coordinates"][0]\n north = coords[1][1]\n south = coords[0][1]\n east = coords[2][0]\n west = coords[0][0]\n center = [(east + west) / 2, (north + south) / 2]\n _meshes[meshcode] = {\n "north": north,\n "south": south,\n "east": east,\n "west": west,\n "center": center\n }\n meshes = _meshes\n\n # 犯罪データの座標から，どのグリッドセルに含まれるかを求め\n # 「mesh_code」列に格納する\n # 読み込んだメッシュGeoJsonに含まれない犯罪データは除外\n df_temp = df_crimes[0:0].copy()\n df_temp["mesh_code"] = None\n for index, row in df_crimes.iterrows():\n for meshcode, coords in meshes.items():\n if ((coords["south"] <= row["latitude"] < coords["north"]) and\n (coords["west"] <= row["longitude"] < coords["east"])):\n if not pd.isnull(row["datetime"]):\n df_temp.loc[index] = pd.concat([\n row,\n pd.Series({\n "mesh_code": meshcode\n })\n ])\n break\n print(f"{df_crimes.index.get_loc(index)} / {len(df_crimes)}")\n df_crimes = df_temp.copy()\n df_crimes["datetime"] = pd.to_datetime(df_crimes["datetime"], format = "%Y/%m/%d %H:%M")\n\n\n # CSVファイルの読み込み\n df_crimes = pd.read_csv(HERE / "犯罪_メッシュ入り.csv")\n #メッシュごとに件数をカウント\n mesh_counts = df_crimes['mesh_code'].value_counts().reset_index()\n mesh_counts.columns = ['mesh_code', '件数']\n # 結果を表示\n print(mesh_counts)\n # 結果をCSVファイルに保存する場合\n mesh_counts.to_csv(HERE / "犯罪_メッシュ件数.csv", index=False, encoding='utf-8-sig')\n\n&br;\n # 件数データ\n df_number = pd.read_csv(\n HERE / "犯罪_メッシュ件数.csv",\n index_col = "mesh_code")\n\n&br;\n # 取引データ\n df_torihiki = pd.read_csv(\n HERE / "df_housing_toyamacityメッシュ入り.csv",\n index_col = "mesh_code",\n encoding = 'cp932')\n\n&br;\n # メッシュコードをキーとして結合\n df_merged = pd.merge(df_torihiki, df_number, on='mesh_code', how='left', \n suffixes=('', '_number'))\n\n\n print(df_merged.head())\n # 結果をCSVファイルに保存する場合\n df_merged.to_csv(HERE / "取引_件数付き.csv", encoding='utf-8-sig')\n \n&br;\n df_mesh = pd.read_csv(\n HERE / "data_by_mesh.csv",\n index_col = "mesh_code",\n encoding = "cp932")\n\n&br;\n df_torihiki2 = pd.read_csv(\n HERE / "取引_件数付き.csv",\n index_col = "mesh_code")\n\n&br;\n # メッシュコードでデータを結合\n merged_df = pd.merge(df_torihiki2, df_mesh, on='mesh_code', how='left')\n print(merged_df.head())\n merged_df.to_csv(HERE / "取引_統計データ付き.csv", index=False, \n encoding='utf-8-sig')\n\n犯罪発生率(件/250000)はhop_メッシュ.ipynbを行った後にcsvファイル上の新しい列(犯罪発生率(件/250000)を追加して行う。\n計算方法は=件数/250000で行う。\n\n#ref(統計データ.png,,20%)\n\n「''distance.ipynb''」を実行すると、住宅地データ (取引_統計データ付き.csv) と犯罪データ (犯罪_メッシュ入り.csv) を読み込み、各住宅地から最寄りの犯罪発生地点までの最短距離を計算するために、geopy.distanceモジュールのgeodesic関数を使用します。この関数は、2点間の地球上での距離（メートル単位）を計算します。契約された土地から最短の距離で犯罪が起きた場所までの最短距離を「''取引_統計データ付き.csv''」に追加した「''住宅地_犯罪_最短距離.csv''」を生成します。\n import pandas as pd\n from geopy.distance import geodesic\n # ファイルを読み込む\n df_housing = pd.read_csv("C:/Users/nt011/Desktop/統計正しいほう/データセット \n の作成/取引_統計データ付き.csv")\n df_crimes = pd.read_csv("C:/Users/nt011/Desktop/統計正しいほう/データセット \n の作成/犯罪_メッシュ入り.csv")\n print(df_housing.columns)\n # 1. 住宅地データの緯度経度情報を抽出\n df_housing['latitude'] = df_housing['緯度']\n df_housing['longitude'] = df_housing['経度']\n # 2. 犯罪データの緯度経度情報を抽出\n df_crimes['latitude'] = df_crimes['latitude']\n df_crimes['longitude'] = df_crimes['longitude']\n # 3.最短距離を計算する関数\n def calculate_min_distance(lat1, lon1, df_target, lat_col, lon_col):\n distances = df_target.apply(lambda row: geodesic((lat1, lon1), \n (row[lat_col], row[lon_col])).meters, axis=1)\n return distances.min()\n # 4. 各住宅地から犯罪発生場所までの最短距離を計算\n df_housing['min_distance_to_crime'] = df_housing.apply(lambda row: \n calculate_min_distance(row['latitude'], row['longitude'], df_crimes, 'latitude', 'longitude'), axis=1)\n # 距離が0になる住宅地の行を表示\n df_housing[df_housing['min_distance_to_crime'] == 0]\n # 必要なカラムを選択して新しいDataFrameに保存\n df_result = df_housing[['mesh_code', '住所', '経度', '緯度', '最寄駅：名称', '最寄駅：距離（分）', '取引価格（総額）', '面積（㎡）', '取引価格（㎡単価）', \n '土地の形状', '間口', '延床面積（㎡）','建築年数', '建物の構造', \n '用途', '前面道路：方位', '前面道路：種類', '前面道路：幅員（ｍ）', '都市計画', '建ぺい率（％）', \n '容積率（％）', '取引時点', '件数', '犯罪発生率(件/250000)', '人口', '18歳未満人口割合', \n '65歳以上人口割合', '外国人人口割合', '世帯数', '単身世帯割合', '核家族世帯割合', '正規労働者割合', \n '非正規労働者割合', '最終学歴が中学以下の人口割合', '最終学歴が高校の人口割合', '最終学歴が大学以上の人口割合', \n '居住年数5年未満の人口割合', '居住年数が20年以上の人口割合', '一戸建て割合', 'アパート・低中層マンション割合', \n '高層マンション割合', 'コンビニ', '駅', '駐車場', '駐輪場', '金融機関', '旅館', 'ホテル', 'スーパー', 'ショッピングモール', 'デパート', '警察署／交番', '小学校', '中学校', '高校', '大学',\n '保育園／幼稚園', 'レジャー施設','道路の面積比率', '建物の面積比率', '水の面積比率', '空き地の面積比率', \n 'ノード数（道路）', 'エッジ数（道路）', '密度（道路）', '平均次数（道路）', 'min_distance_to_crime']]\n # 結果の保存\n df_result.to_csv('C:/Users/nt011/Desktop/統計正しいほう/データセットの作成/ \n 住宅地_犯罪_最短距離.csv', index=False, encoding='utf-8-sig')\n\n#ref(最短.png,,20%)&br;\n&br;\n** データセットの前処理 [#u50edad0]\n\n作成したデータセットを回帰分析がかけられる形に整える。\n&br;\n\n「''前処理.ipynb''」を実行し、「''住宅地_犯罪_最短距離.csv''」から「建築年数」と「最寄駅：距離（分）」の列に関する欠損値を含む行を削除した「''住宅地_犯罪_最短距離_欠損値削除後.csv''」を生成します。\n&br;\n import pandas as pd\n # CSVファイルのパス\n file_path = r"C:\Users\nt011\Desktop\統計正しいほう\前処理\住宅地_犯罪_最短距離.csv"\n # CSVファイルを読み込む\n df = pd.read_csv(file_path)\n # 欠損値を含む行を削除\n df = df.dropna()\n # 「建築年数」列に `-1` または `"VALUE"` を含む行を削除\n df = df[~df["建築年数"].astype(str).isin(["-1", "#VALUE!"])]\n # 「最寄駅：距離（分）」の値を変換\n replace_dict = {\n "30分~60分": 45,\n "1H~1H30": 75,\n "1H30~2H": 105\n }\n df["最寄駅：距離（分）"] = df["最寄駅：距離（分）"].replace(replace_dict)\n # クリーンなデータを新しいCSVファイルとして保存\n output_path = r"C:\Users\nt011\Desktop\統計正しいほう\前処理\住宅地_犯罪_最短距離_欠損値削除後.csv"\n df.to_csv(output_path, index=False, encoding='utf-8-sig')\n print("データ処理完了: 欠損値と指定された異常値を削除し、駅距離の値を変換しました。")\n print(f"新しいCSVファイルが {output_path} に保存されました。")\n\n\n生成した「''住宅地_犯罪_最短距離_欠損値削除後.csv''」から連続変数の外れ値を検出するカラムリストを定義し、データ型を数値に変換したのち検出された列を削除したリスト「''外れ値削除.csv''」を生成します。ダミー変数の処理を行う前にcsvファイル内の建物の構造の英語の部分を変換する必要がある。\n\n◆不動産取引価格~\nhttps://www.reinfolib.mlit.go.jp/realEstatePrices/~\n~\nこのurlの制度と用語のページを参考にして変換を行う。\n#ref(建物の構造.png,,30%)\n&br;\n import pandas as pd\n # CSVファイルのパス\n file_path = "C:/Users/nt011/Desktop/統計正しいほう/前処理/住宅地_犯罪_最短距 \n 離_欠損値削除後.csv"\n # CSVファイルを読み込む\n df = pd.read_csv(file_path)\n # 外れ値を検出する連続変数のカラムリスト\n continuous_columns = [\n "最寄駅：距離（分）", "取引価格（総額）", "面積（㎡）", "取引価格（㎡単価）",\n "間口", "延床面積（㎡）", "建築年数", "前面道路：幅員（ｍ）",\n "建ぺい率（％）", "容積率（％）", "犯罪発生率(件/250000)", "人口",\n "18歳未満人口割合", "65歳以上人口割合", "外国人人口割合", "世帯数",\n "単身世帯割合", "核家族世帯割合", "正規労働者割合", "非正規労働者割合",\n "最終学歴が中学以下の人口割合", "最終学歴が高校の人口割合", "最終学歴が大学以上の人口割合",\n "居住年数5年未満の人口割合", "居住年数が20年以上の人口割合", "一戸建て割合",\n "アパート・低中層マンション割合", "高層マンション割合", "道路の面積比率",\n "建物の面積比率", "水の面積比率", "空き地の面積比率", "ノード数（道路）",\n "エッジ数（道路）", "密度（道路）", "平均次数（道路）", "min_distance_to_crime"\n ]\n # データ型を数値に変換（エラーが出たらNaNにする）\n for column in continuous_columns:\n if column in df.columns:\n df[column] = pd.to_numeric(df[column], errors='coerce')\n # 外れ値を持つ行のインデックスを格納するリスト\n outlier_indices = set()\n # 外れ値の検出と削除\n for column in continuous_columns:\n if column in df.columns:\n # 欠損値を除外\n df_clean = df[column].dropna()\n # IQR（四分位範囲）の計算\n Q1 = df_clean.quantile(0.25)\n Q3 = df_clean.quantile(0.75)\n IQR = Q3 - Q1\n # 外れ値の閾値を設定\n lower_bound = Q1 - 1.5 * IQR\n upper_bound = Q3 + 1.5 * IQR\n # 外れ値のインデックスを取得\n outliers = df[(df[column] < lower_bound) | (df[column] > upper_bound)]\n outlier_indices.update(outliers.index)\n # 外れ値の行を削除\n df_cleaned = df.drop(index=outlier_indices)\n # 結果を新しいCSVに保存\n output_path = "C:/Users/nt011/Desktop/統計正しいほう/前処理/外れ値削除.csv"\n df_cleaned.to_csv(output_path, encoding="utf-8-sig", index=False)\n print(f"外れ値を削除したデータを {output_path} に保存しました。")\n print(f"削除した行数: {len(outlier_indices)}")\n print(f"最終データの行数: {df_cleaned.shape[0]}")\n\n\nそれぞれのカラムの要素を必要に応じてダミー変数化した後、データフレームを保存した「''削除_ダミー後.csv''」を生成します。~\n import pandas as pd\n # CSVファイルを読み込み\n df_data = pd.read_csv("C:/Users/nt011/Desktop/統計正しいほう/前処理/外れ値削除.csv")\n # カラム名を取得\n columns = df_data.columns.tolist()\n # 建物の構造列に「RC」がある場合、「鉄筋コンクリート造」に変換\n df_data['建物の構造'] = df_data['建物の構造'].replace('RC', '鉄筋コンクリート造')\n # 最寄り駅名を二値化（列名を「最寄り駅名_富山」に変更）\n index = columns.index("最寄駅：名称") # 変換前の位置を取得\n df_data.insert(index + 1, '最寄り駅名_富山', df_data['最寄駅：名称'].apply(lambda x: 1 if x == '富山' else 0))\n df_data = df_data.drop(columns=['最寄駅：名称']) # 元の列を削除\n # 都市計画を二値化（列名を「都市計画_市街化調整区域」に変更）\n index = columns.index("都市計画")\n df_data.insert(index + 1, '都市計画_市街化調整区域', df_data['都市計画'].apply(lambda x: 1 if x == '市街化調整区域' else 0))\n df_data = df_data.drop(columns=['都市計画'])\n # 土地の形状を二値化（列名を「土地の形状_不整形」に変更）\n index = columns.index("土地の形状")\n df_data.insert(index + 1, '土地の形状_不整形', df_data['土地の形状'].apply(lambda x: 1 if x == '不整形' else 0))\n df_data = df_data.drop(columns=['土地の形状'])\n # 用途を二値化（列名を「用途_住宅」に変更）\n index = columns.index("用途")\n df_data.insert(index + 1, '用途_住宅', df_data['用途'].apply(lambda x: 1 if x == '住宅' else 0))\n df_data = df_data.drop(columns=['用途'])\n # 建物の構造を二値化（列名を「建物の構造_木造」に変更）\n index = columns.index("建物の構造")\n df_data.insert(index + 1, '建物の構造_木造', df_data['建物の構造'].apply(lambda x: 1 if x == '木造' else 0))\n df_data = df_data.drop(columns=['建物の構造'])\n # 「前面道路：方位」のダミー変数化（基準値「北西」を0に設定）\n index = columns.index("前面道路：方位") # 変換前の位置を取得\n dummies = pd.get_dummies(df_data['前面道路：方位'], prefix='前面道路：方位', dtype=int)\n # 基準値は削除\n baseline_value = '北西'\n if f'前面道路：方位_{baseline_value}' in dummies.columns:\n dummies = dummies.drop(columns=[f'前面道路：方位_{baseline_value}'])\n # 挿入（前面道路：方位のダミー変数）\n ordered_dummies = ['前面道路：方位_西', '前面道路：方位_東', '前面道路：方位_南西', '前面道路：方位_南東', \n '前面道路：方位_南', '前面道路：方位_北東', '前面道路：方位_北']\n # 順番に従って挿入\n for col in ordered_dummies:\n if col in dummies.columns:\n df_data.insert(index + 1, col, dummies[col])\n # 「前面道路：種類」のダミー変数化\n road_types = ['市道', '私道', '県道']\n # 順番を指定してダミー変数を挿入\n for road in road_types:\n df_data.insert(index + 1, f'前面道路：種類_{road}', (df_data['前面道路：種類'] == road).astype(int))\n # 元の「前面道路：方位」および「前面道路：種類」列を削除\n df_data = df_data.drop(columns=['前面道路：方位', '前面道路：種類'])\n # 確認用にダミー変数化後のデータフレームを表示\n print(df_data.head())\n # ダミー変数化後のデータを保存\n df_data.to_csv("C:/Users/nt011/Desktop/統計正しいほう/前処理/削除_ダミー後.csv", \n index=False, encoding='utf-8-sig')\n print("ダミー変数変換完了！データが保存されました。")\n\n#ref(ダミー.png,,20%)\n&br;\n\n連続変数とダミー変数を手動で分類し、連続変数のみ二次項の作成をする。\n交互作用項の作成したデータフレームを保存した「''削除_変数作成後.csv''」を生成します。\n import pandas as pd\n import numpy as np\n from itertools import combinations\n # CSVファイルのパスを指定\n file_path = "C:/Users/nt011/Desktop/統計正しいほう/前処理/削除_ダミー後.csv"\n # UTF-8で読み込む\n df = pd.read_csv(file_path)\n # データの先頭を表示\n print(df.head())\n # 連続変数とダミー変数を手動で分類\n continuous_vars = [\n '最寄駅：距離（分）', '面積（㎡）', '間口', '延床面積（㎡）', '建築年数', '前面道路：幅員（ｍ）',\n '建ぺい率（％）', '容積率（％）', '犯罪発生率(件/250000)', '人口', '18歳未満人口割合',\n '65歳以上人口割合', '外国人人口割合', '世帯数', '単身世帯割合', '核家族世帯割合', \n '正規労働者割合', '非正規労働者割合', '最終学歴が中学以下の人口割合', '最終学歴が高校の人口割合',\n '最終学歴が大学以上の人口割合', '居住年数5年未満の人口割合', '居住年数が20年以上の人口割合',\n '一戸建て割合', 'アパート・低中層マンション割合', 'コンビニ', '駅', '駐車場', '駐輪場', \n '金融機関', '旅館', 'ホテル', 'スーパー', 'ショッピングモール', \n 'デパート', '警察署／交番', '小学校', '中学校', '高校', '大学', '保育園／幼稚園', \n 'レジャー施設','道路の面積比率', '建物の面積比率', '水の面積比率', '空き地の面積比率', \n 'ノード数（道路）', 'エッジ数（道路）', '密度（道路）', '平均次数（道路）', 'min_distance_to_crime'\n ]\n # ダミー変数をリストとして定義\n dummy_vars = [\n '最寄り駅名_富山','土地の形状_不整形', '建物の構造_木造', \n '用途_住宅', '前面道路：方位_北', '前面道路：方位_北東', \n '前面道路：方位_南', '前面道路：方位_南東', \n '前面道路：方位_南西', '前面道路：方位_東', '前面道路：方位_西', \n '前面道路：種類_市道', '前面道路：種類_県道', '前面道路：種類_私道', \n '都市計画_市街化調整区域'\n ]\n # 二次項の作成（連続変数のみ）\n for var in continuous_vars:\n df[f'{var}^2'] = df[var] ** 2\n # 交互作用項の作成\n # 連続変数 × 連続変数\n for var1, var2 in combinations(continuous_vars, 2):\n df[f'{var1}×{var2}'] = df[var1] * df[var2]\n # ダミー変数 × 連続変数\n for dummy_var in dummy_vars:\n for cont_var in continuous_vars:\n df[f'{dummy_var}×{cont_var}'] = df[dummy_var] * df[cont_var]\n # ダミー変数 × ダミー変数\n for var1, var2 in combinations(dummy_vars, 2):\n df[f'{var1}×{var2}'] = df[var1] * df[var2]\n # CSVに保存\n output_path = "C:/Users/nt011/Desktop/統計正しいほう/前処理/削除_変数作成後.csv"\n df.to_csv(output_path, index=False, encoding='utf-8-sig')\n print(f"処理が完了しました。生成されたデータは {output_path} に保存されました。")\n#ref(交互作用項.png,,20%)\n&br;\n\n*変数選択 [#u9327c42]\n***Elastic Netとは [#d2ec2870]\nElastic Netは、L1 正則化（Lasso）と L2 正則化（Ridge）の組み合わせにより、多重共線性があるデータに対して安定した回帰分析を行う手法です。\n\nLasso（L1正則化）：L1正則化は、回帰モデルの係数の一部をゼロにすることによって、特徴量選択を行います。これにより、重要な特徴量だけを残すことができますが、多重共線性が高い場合には不安定になることがあります。\n&br;\n\nRidge（L2正則化）：L2正則化は、回帰モデルの係数を小さくすることで、過学習を防ぎます。L2正則化は係数をゼロにすることはなく、むしろ係数の大きさを均等に抑える傾向がありますが、特徴量選択の機能はありません。\n&br;\n\nStandardScaler() を用いて X（説明変数）と y（目的変数）を標準化 しているため、変数のスケールによる影響を排除し、モデルの安定性を向上させています。\n\n\nAdaptive Elastic Netを用いて土地の価格に影響する重要な特徴量の係数を取得し、影響の大きい特徴量を特定します。&br;\n#ref(式.png,,50%)\n***手順 [#d2ec2870]\n「''変数選択__削除.ipynb''」を実行すると、重要な変数のみを抽出し、「''選ばれた変数.csv''」に保存します。\n import numbers\n import warnings\n import cvxpy\n import numpy as np\n from asgl import ASGL\n from sklearn.base import MultiOutputMixin\n from sklearn.base import RegressorMixin\n from sklearn.exceptions import ConvergenceWarning\n from sklearn.linear_model import ElasticNet\n from sklearn.linear_model._coordinate_descent import _alpha_grid\n from sklearn.utils import check_X_y\n from sklearn.utils.validation import check_is_fitted\n class AdaptiveElasticNet(ASGL, ElasticNet, MultiOutputMixin, RegressorMixin):\n """\n Objective function and parameters as described with modifications\n to allow alpha1, alpha2, l1_ratio1, and l1_ratio2 customization.\n """\n def __init__(\n self,\n alpha1=0.021544346900318846, # First-stage ElasticNet alpha\n alpha2=0.0009443498043343188, # Second-stage ElasticNet alpha\n *,\n l1_ratio1=0.8889, # First-stage ElasticNet L1 ratio\n l1_ratio2=0.778, # Second-stage ElasticNet L1 ratio\n gamma=0.5, # Weight adjustment exponent\n fit_intercept=True,\n precompute=False,\n max_iter=10000,\n copy_X=True,\n solver=None,\n tol=None,\n positive=False,\n positive_tol=1e-3,\n random_state=None,\n eps_coef=1e-6,\n verbose=True\n ):\n params_asgl = dict(model="lm", penalization="asgl")\n if solver is not None:\n params_asgl["solver"] = solver\n if tol is not None:\n params_asgl["tol"] = tol\n super().__init__(**params_asgl)\n self.alpha1 = alpha1\n self.alpha2 = alpha2\n self.l1_ratio1 = l1_ratio1\n self.l1_ratio2 = l1_ratio2\n self.gamma = gamma\n self.fit_intercept = fit_intercept\n self.max_iter = max_iter\n self.precompute = precompute\n self.copy_X = copy_X\n self.positive = positive\n self.positive_tol = positive_tol\n self.random_state = random_state\n self.eps_coef = eps_coef\n self.verbose = verbose\n if not self.fit_intercept:\n raise NotImplementedError\n def fit(self, X, y, check_input=True):\n if check_input:\n X_copied = self.copy_X and self.fit_intercept\n X, y = self._validate_data(\n X,\n y,\n accept_sparse="csc",\n order="F",\n dtype=[np.float64, np.float32],\n copy=X_copied,\n multi_output=True,\n y_numeric=True,\n )\n # 第一段階の ElasticNet 実行\n enet_model_1 = self.elastic_net(self.l1_ratio1, alpha=self.alpha1)\n enet_model_1.fit(X, y)\n enet_coef = enet_model_1.coef_\n # 重みの計算\n weights = 1.0 / (np.maximum(np.abs(enet_coef), self.eps_coef) ** self.gamma)\n # 第二段階の最適化\n self.coef_, self.intercept_ = self._optimize_second_stage(X, y, weights)\n # モデル属性を格納\n self.enet_coef_ = enet_coef\n self.weights_ = weights\n return self\n def predict(self, X):\n check_is_fitted(self, ["coef_", "intercept_"])\n return super(ElasticNet, self).predict(X)\n def elastic_net(self, l1_ratio, **params):\n """\n Create an ElasticNet model with the specified parameters.\n """\n elastic_net = ElasticNet(l1_ratio=l1_ratio)\n for k, v in self.get_params().items():\n try:\n elastic_net = elastic_net.set_params(**{k: v})\n except ValueError:\n pass # Ignore parameters not supported by ElasticNet\n elastic_net = elastic_net.set_params(**params)\n return elastic_net\n def _optimize_second_stage(self, X, y, weights):\n """\n Perform second-stage optimization with adaptive weights.\n Returns\n -------\n coef : np.array, shape (n_features,)\n intercept : float\n """\n n_samples, n_features = X.shape\n beta_variables = [cvxpy.Variable(n_features)]\n model_prediction = 0.0\n if self.fit_intercept:\n beta_variables = [cvxpy.Variable(1)] + beta_variables\n ones = cvxpy.Constant(np.ones((n_samples, 1)))\n model_prediction += ones @ beta_variables[0]\n # モデル予測\n model_prediction += X @ beta_variables[1]\n error = cvxpy.sum_squares(y - model_prediction) / (2 * n_samples)\n l1_coefs = self.alpha2 * self.l1_ratio2\n # 第二段階の正則化項\n l1_penalty = cvxpy.Constant(l1_coefs * weights) @ cvxpy.abs(\n beta_variables[1]\n )\n l2_penalty = (\n cvxpy.Constant(self.alpha1 * (1 - self.l1_ratio1))\n * cvxpy.sum_squares(beta_variables[1])\n )\n constraints = [b >= 0 for b in beta_variables] if self.positive else []\n # 最適化問題の定義\n problem = cvxpy.Problem(\n cvxpy.Minimize(error + l1_penalty + l2_penalty), \n constraints=constraints\n )\n problem.solve(solver="OSQP", max_iter=self.max_iter)\n if problem.status != "optimal":\n raise ConvergenceWarning(\n f"Solver did not reach optimum (Status: {problem.status})"\n )\n beta_sol = np.concatenate([b.value for b in beta_variables], axis=0)\n beta_sol[np.abs(beta_sol) < self.tol] = 0\n intercept, coef = beta_sol[0], beta_sol[1:]\n coef = np.maximum(coef, 0) if self.positive else coef\n return coef, intercept\n from sklearn.linear_model import ElasticNetCV\n from sklearn.model_selection import GridSearchCV\n from sklearn.datasets import make_regression\n import numpy as np\n from sklearn.model_selection import train_test_split\n import pandas as pd\n from sklearn.preprocessing import StandardScaler\n # CSVファイルのパスを指定してください\n file_path = "C:/Users/nt011/Desktop/統計正しいほう/変数選択/削除_変数作成後.csv"\n # CSVファイルの読み込み\n df = pd.read_csv(file_path) # 文字コードが異なる場合は、'utf-8' を他のエンコーディングに変更してください\n # 特徴量とターゲットの分割\n X = df.drop(columns=['取引価格（㎡単価）']) # 取引価格（㎡単価）をyとして分離\n y = df['取引価格（㎡単価）'] # ターゲット変数\n # 🔥 標準化の実施\n scaler_X = StandardScaler()\n scaler_y = StandardScaler()\n X_scaled = scaler_X.fit_transform(X) # Xの標準化\n y_scaled = scaler_y.fit_transform(y.values.reshape(-1, 1)).ravel() # yの標準化\n # データ分割\n X_train, X_test, y_train, y_test = train_test_split(X_scaled, y_scaled, \n test_size=0.2, random_state=42)\n # 訓練データの一部を使用\n X_sample, _, y_sample, _ = train_test_split(X_train, y_train, test_size=0.8, \n random_state=42)\n # 第一段階: ElasticNetCVを使用した最適パラメータの導出\n enet_cv = ElasticNetCV(\n l1_ratio=np.linspace(0.0001, 1, 25), # l1_ratioの候補\n alphas=np.logspace(-5, 0, 25), # alphaの候補\n cv=5, # 交差検証の分割数\n random_state=42,\n n_jobs=-1,\n verbose=1\n )\n enet_cv.fit(X_train, y_train)\n # 第一段階の最適パラメータと係数を取得\n alpha1_opt = enet_cv.alpha_\n l1_ratio1_opt = enet_cv.l1_ratio_\n enet_coef = enet_cv.coef_\n print(f"第一段階の最適パラメータ: alpha1={alpha1_opt}, l1_ratio1= \n {l1_ratio1_opt}")\n\n # ベストパラメータをalpha2とl1_ratio2でそれぞれ格納\n best_alpha2 = best_params['alpha2']\n best_l1_ratio2 = best_params['l1_ratio2']\n print(best_alpha2)\n print(best_l1_ratio2)\nASGLとElasticNetを組み合わせたAdaptiveElasticNetクラスを定義し、CSVから読み込んだデータを標準化後にモデルをフィッティングして、非ゼロの係数を持つ重要な特徴量を抽出し、その結果をCSVに保存する一連の処理を実行します。\n import numbers\n import warnings\n import cvxpy\n import numpy as np\n from asgl import ASGL\n from sklearn.base import MultiOutputMixin\n from sklearn.base import RegressorMixin\n from sklearn.exceptions import ConvergenceWarning\n from sklearn.linear_model import ElasticNet\n from sklearn.linear_model._coordinate_descent import _alpha_grid\n from sklearn.utils import check_X_y\n from sklearn.utils.validation import check_is_fitted\n class AdaptiveElasticNet(ASGL, ElasticNet, MultiOutputMixin, RegressorMixin):\n """\n Objective function and parameters as described with modifications\n to allow alpha1, alpha2, l1_ratio1, and l1_ratio2 customization.\n """\n def __init__(\n self,\n alpha1=0.01333521432163324, # 第一段階で得られたalpha1を指定\n alpha2=0.034807005884284134, # Second-stage ElasticNet alpha\n *,\n l1_ratio1=1, # First-stage ElasticNet L1 ratio\n l1_ratio2=0.1667499999999998, # Second-stage ElasticNet L1 ratio\n gamma=0.5, # Weight adjustment exponent\n fit_intercept=True,\n precompute=False,\n max_iter=10000,\n copy_X=True,\n solver=None,\n tol=None,\n positive=False,\n positive_tol=1e-3,\n random_state=None,\n eps_coef=1e-6,\n verbose=True\n ):\n params_asgl = dict(model="lm", penalization="asgl")\n if solver is not None:\n params_asgl["solver"] = solver\n if tol is not None:\n params_asgl["tol"] = tol\n super().__init__(**params_asgl)\n self.alpha1 = alpha1\n self.alpha2 = alpha2\n self.l1_ratio1 = l1_ratio1\n self.l1_ratio2 = l1_ratio2\n self.gamma = gamma\n self.fit_intercept = fit_intercept\n self.max_iter = max_iter\n self.precompute = precompute\n self.copy_X = copy_X\n self.positive = positive\n self.positive_tol = positive_tol\n self.random_state = random_state\n self.eps_coef = eps_coef\n self.verbose = verbose\n if not self.fit_intercept:\n raise NotImplementedError\n def fit(self, X, y, check_input=True):\n if check_input:\n X_copied = self.copy_X and self.fit_intercept\n X, y = self._validate_data(\n X,\n y,\n accept_sparse="csc",\n order="F",\n dtype=[np.float64, np.float32],\n copy=X_copied,\n multi_output=True,\n y_numeric=True,\n )\n # 第一段階の ElasticNet 実行\n enet_model_1 = self.elastic_net(self.l1_ratio1, alpha=self.alpha1)\n enet_model_1.fit(X, y)\n enet_coef = enet_model_1.coef_\n # 重みの計算\n weights = 1.0 / (np.maximum(np.abs(enet_coef), self.eps_coef) ** self.gamma)\n # 第二段階の最適化\n self.coef_, self.intercept_ = self._optimize_second_stage(X, y, weights)\n # モデル属性を格納\n self.enet_coef_ = enet_coef\n self.weights_ = weights\n return self\n def predict(self, X):\n check_is_fitted(self, ["coef_", "intercept_"])\n return super(ElasticNet, self).predict(X)\n def elastic_net(self, l1_ratio, **params):\n """\n Create an ElasticNet model with the specified parameters.\n """\n elastic_net = ElasticNet(l1_ratio=l1_ratio)\n for k, v in self.get_params().items():\n try:\n elastic_net = elastic_net.set_params(**{k: v})\n except ValueError:\n pass # Ignore parameters not supported by ElasticNet\n elastic_net = elastic_net.set_params(**params)\n return elastic_net\n def _optimize_second_stage(self, X, y, weights):\n """\n Perform second-stage optimization with adaptive weights.\n Returns\n -------\n coef : np.array, shape (n_features,)\n intercept : float\n """\n n_samples, n_features = X.shape\n beta_variables = [cvxpy.Variable(n_features)]\n model_prediction = 0.0\n if self.fit_intercept:\n beta_variables = [cvxpy.Variable(1)] + beta_variables\n ones = cvxpy.Constant(np.ones((n_samples, 1)))\n model_prediction += ones @ beta_variables[0]\n # モデル予測\n model_prediction += X @ beta_variables[1]\n error = cvxpy.sum_squares(y - model_prediction) / (2 * n_samples)\n l1_coefs = self.alpha2 * self.l1_ratio2\n # 第二段階の正則化項\n l1_penalty = cvxpy.Constant(l1_coefs * weights) @ cvxpy.abs(\n beta_variables[1]\n )\n l2_penalty = (\n cvxpy.Constant(self.alpha1 * (1 - self.l1_ratio1))\n * cvxpy.sum_squares(beta_variables[1])\n )\n constraints = [b >= 0 for b in beta_variables] if self.positive else []\n # 最適化問題の定義\n problem = cvxpy.Problem(\n cvxpy.Minimize(error + l1_penalty + l2_penalty), \n constraints=constraints\n )\n problem.solve(solver="OSQP", max_iter=self.max_iter)\n if problem.status != "optimal":\n raise ConvergenceWarning(\n f"Solver did not reach optimum (Status: {problem.status})"\n )\n beta_sol = np.concatenate([b.value for b in beta_variables], axis=0)\n beta_sol[np.abs(beta_sol) < self.tol] = 0\n intercept, coef = beta_sol[0], beta_sol[1:]\n coef = np.maximum(coef, 0) if self.positive else coef\n return coef, intercept\n import pandas as pd\n from sklearn.preprocessing import StandardScaler\n # CSVファイルのパスを指定してください\n file_path = "C:/Users/monda/OneDrive/デスクトップ/統計正しいほう/変数選択/削除_変数作成後.csv"\n # CSVファイルの読み込み\n df = pd.read_csv(file_path) # 文字コードが異なる場合は、'utf-8' を他のエンコー ディングに変更してください\n # 特徴量とターゲットの分割\n X = df.drop(columns=['取引価格（㎡単価）']) # 取引価格（㎡単価）をyとして分離\n y = df['取引価格（㎡単価）'] # ターゲット変数\n # 🔥 標準化の実施\n scaler_X = StandardScaler()\n scaler_y = StandardScaler()\n X_scaled = scaler_X.fit_transform(X) # Xの標準化\n y_scaled = scaler_y.fit_transform(y.values.reshape(-1, 1)).ravel() # yの標準化\n model = AdaptiveElasticNet()\n model.fit(X_scaled, y_scaled)\n # 係数の表示（重要な特徴量を確認）\n coef = model.coef_\n print("特徴量の係数:", coef)\n # 係数がゼロでない特徴量を選定\n selected_features = X.columns[coef != 0]\n print("選定された特徴量:", selected_features)\n # 係数と特徴量名をデータフレームに変換\n selected_features_df = pd.DataFrame({\n 'Feature': X.columns,\n 'Coefficient': coef\n })\n # 係数がゼロでない特徴量のみをフィルタリング\n selected_features_df = selected_features_df[selected_features_df['Coefficient'] != 0]\n # CSVファイルとして保存\n output_file_path = "選ばれた変数.csv" # 保存先を指定してください\n selected_features_df.to_csv(output_file_path, index=False, encoding='utf-8-sig') \n # 日本語対応のエンコーディング\n print(f"選定された特徴量と係数を {output_file_path} に保存しました。")\n #ref(変数選択.png,,50%)\n&br;\nAENの数式は以下のようなものを扱う。そして、今回行ったハイパーパラメータは以下のものである。\nAENの推計は，二段階で行われる．まず，予備推計としてEN で係数を推定する．その\nうえで，係数の絶対値が小さい変数に対してより強い罰則を与えるようL1 ノルムの正則化\n項を変数ごとに調整したうえで，改めてEN を実施する．\n#ref(AEN.png,,150%)\n#ref(ハイパーパラメータ.png,,150%)\n\n選択された変数は以下のようになる\n#ref(選ばれた変数.png,,)\n\n*要因分析 [#w249dca6]\n***手順 [#d20580fd]\n「''hedo.ipynb''」を実行すると、有意水準に基づいて有意性を判定する関数\n import pandas as pd\n # Excelファイルを読み込み\n file_path = "c:/Users/nt011/Desktop/統計正しいほう/要因分析/位置情報_削除.csv"\n df = pd.read_csv(file_path)\n # 住所が含まれている列を特定して、住所ごとの出現回数を確認\n address_column = '住所' # 住所列の名前がわからない場合はここを変更してください\n if address_column in df.columns:\n address_counts = df[address_column].dropna().value_counts() # 住所ごとの出現回数を集計\n # 出現回数をデータフレームに変換\n address_counts_df = address_counts.reset_index()\n address_counts_df.columns = ['住所', '出現回数'] # 列名を設定\n # 結果をCSVファイルに書き出し\n output_file = "c:/Users/nt011/Desktop/統計正しいほう/要因分析/住所出現回数_削除.csv"\n address_counts_df.to_csv(output_file, index=False, encoding='utf-8-sig') \n # UTF-8エンコーディングで保存\n print(f"住所ごとの出現回数がCSVファイルに保存されました: {output_file}")\n else:\n print("住所列が見つかりません")\n回帰係数、t値、p値、有意性を格納したリストを「''回帰分析結果_削除.csv''」に生成します。\n import numpy as np\n import pandas as pd\n from sklearn.preprocessing import StandardScaler\n # ファイルパスの指定\n selected_features_path = "選ばれた変数_削除.csv"\n all_data_path = "c:/Users/nt011/Desktop/統計正しいほう/要因分析/削除_変数作成後.csv"\n location_data_path = "c:/Users/nt011/Desktop/統計正しいほう/要因分析/位置情報_削除.csv"\n output_data_path = "選択された特徴量とデータ_削除.csv"\n # 選定された変数のCSVファイルを読み込み\n selected_features_df = pd.read_csv(selected_features_path, encoding='utf-8-sig')\n selected_features = selected_features_df['Feature'].tolist()\n # 元データを再度読み込み\n all_data_df = pd.read_csv(all_data_path)\n # 選定された説明変数と目的変数を抽出\n df = all_data_df[selected_features + ['取引価格（㎡単価）']]\n # 位置情報データを読み込み\n location_df = pd.read_csv(location_data_path, usecols=['住所'])\n # 元データと位置情報を結合（共通のキーがある場合）\n df = pd.concat([df, location_df], axis=1)\n # 統合データをCSVファイルとして保存\n df.to_csv(output_data_path, index=False, encoding='utf-8-sig')\n print(f"選定された変数と位置情報を含むデータを {output_data_path} に保存しました。")\n\n import pandas as pd\n # ファイルパス\n count_file = "c:/Users/nt011/Desktop/統計正しいほう/要因分析/住所出現回数_削除.csv"\n data_file = "C:/Users/nt011/Desktop/統計正しいほう/要因分析/選択された特徴量とデータ_削除.csv"\n output_file = "C:/Users/nt011/Desktop/統計正しいほう/要因分析/ソート済みデータ_削除.csv"\n # CSVファイルの読み込み\n count_df = pd.read_csv(count_file, encoding='utf-8-sig')\n data_df = pd.read_csv(data_file, encoding='utf-8-sig')\n # 住所列の名前\n address_column = '住所'\n # 緯度・経度を除外（該当カラムがある場合）\n drop_columns = ['緯度', '経度']\n data_df = data_df.drop(columns=[col for col in drop_columns if col in \n data_df.columns])\n # 住所の出現回数をデータに統合\n merged_df = data_df.merge(count_df, on=address_column, how='left')\n # 出現回数の多い順にソート\n sorted_df = merged_df.sort_values(by='出現回数', ascending=False).drop(columns=['出現回数'])\n # ソート済みデータをCSVに保存\n sorted_df.to_csv(output_file, encoding='utf-8-sig', index=False)\n print(f"ソート済みデータを {output_file} に保存しました。")\n import pandas as pd\n import numpy as np\n from sklearn.linear_model import LinearRegression\n from sklearn.preprocessing import StandardScaler\n import statsmodels.api as sm\n # データの読み込み\n file_path = "C:/Users/nt011/Desktop/統計正しいほう/要因分析/ソート済みデータ_削除.csv"\n data = pd.read_csv(file_path)\n # 目的変数と説明変数の設定\n y = data['取引価格（㎡単価）'] # 取引価格を目的変数に\n X = data.drop(columns=['取引価格（㎡単価）', '住所']) # 取引価格と住所以外を説明変数に\n # 🔥 標準化の実施\n scaler_X = StandardScaler()\n scaler_y = StandardScaler()\n X_scaled = scaler_X.fit_transform(X) # Xの標準化\n y_scaled = scaler_y.fit_transform(y.values.reshape(-1, 1)).ravel() # yの標準化\n # 線形回帰モデルの作成\n model = LinearRegression()\n model.fit(X_scaled, y)\n # 回帰分析の結果を表示\n X_scaled_with_const = sm.add_constant(X_scaled) # 定数項を追加\n ols_model = sm.OLS(y, X_scaled_with_const).fit()\n # 必要な情報を抽出\n r_squared = ols_model.rsquared # 決定係数\n coefficients = ols_model.params # 回帰係数（定数項も含む）\n t_values = ols_model.tvalues # t値（定数項も含む）\n p_values = ols_model.pvalues # p値（有意水準を判定するため）\n variable_names = ['定数項'] + list(X.columns)\n # 有意水準に基づいて有意性を判定する関数\n def significance(p_value):\n if p_value < 0.01:\n return "1% 有意"\n elif p_value < 0.05:\n return "5% 有意"\n elif p_value < 0.10:\n return "10% 有意"\n else:\n return "有意でない"\n # 結果を格納するリスト\n results = []\n # 回帰係数、t値、p値、有意性をリストに格納\n for var_name, coef, t_val, p_val in zip(variable_names, coefficients, t_values, p_values):\n significance_level = significance(p_val)\n results.append([var_name, coef, t_val, p_val, significance_level])\n # 結果をDataFrameに変換\n results_df = pd.DataFrame(results, columns=["変数", "回帰係数", "t値", "p値", "有意性"])\n # 予測値の計算\n data["標準化予測取引価格（㎡単価）"] = ols_model.predict(X_scaled_with_const)\n data["標準化実測取引価格（㎡単価）"] = y_scaled\n # **結果を表示**\n print("\n==== 回帰分析結果 ====")\n print(f"決定係数 (R^2): {r_squared:.4f}\n")\n print(results_df)\n # 予測結果を含めたデータを保存\n predicted_output_file_path = "C:/Users/nt011/Desktop/統計正しいほう/要因分析/予測結果_削除.csv"\n data.to_csv(predicted_output_file_path, index=False, encoding="utf-8-sig")\n # CSVファイルとして保存\n output_file_path = "C:/Users/nt011/Desktop/統計正しいほう/要因分析/回帰分析結果_削除.csv"\n results_df.to_csv(output_file_path, index=False, encoding="utf-8-sig")\n print(f"回帰分析結果が{output_file_path}に保存されました。")\n#ref(ヘドニック.png,,50%)\n)
タイムスタンプを変更しない
[[中島]] ---- ~ &size(15){''[[参考文献とかメモとか>./参考文献とかメモとか]]''};~ ~ ---- #memo([[技術資料]]\n*目次 [#f3c34276]\n\n#CONTENTS\n\n*本研究の目的 [#cec342b3]\n本研究の目的は，ヘドニック・アプローチにおける限界を克服し，価格形成におけるより\n正確で信頼性の高い評価手法を提供することである．ヘドニックモデルでは，多重共線性，\n欠落変数のバイアス，交互作用効果および非線形効果の取り扱いに問題がある[2]．これら\nは不動産市場や消費財の価格評価において誤差を生じさせる要因となる．本研究は，これ\nらの問題を解決するために，スパース推定法を活用した新たなアプローチを提案する．\n従来のヘドニックモデルは，そのシンプルな線形構造と主要変数の関係を明確に示す点\nで強力であったが，現実の市場では非線形な相互作用や高次の効果が絡むため，これらを\n十分に捉えることができないことが多かった．また，異なる要因が相互作用を持つ場合，そ\nの効果を正確に評価するためには，より複雑なモデルが求められる．そこで本研究では，\nElasticNetおよびその進化形であるAdaptiveElasticNet(AEN)を用い，従来の限界を克服\nすることを目指す．\nElasticNetはL1およびL2正則化を組み合わせることで，多重共線性の問題に対応する\nと同時に，重要な変数の選択を行う能力を持つ．これにより，データ中で重要な特徴を捉\nえ，効果的な変数選択を実現できる．また，AENは，従来のElasticNetをさらに発展させ，\nデータの特性に応じて正則化パラメータを動的に調整する．この適応的なアプローチは，複\n雑なデータ構造に対して非常に効果的であり，非線形性や交互作用を持つデータにも適用\nできる．これにより，不動産市場や消費財の価格形成における微細な変動や相互作用を捉\nえ，価格予測の精度を向上させることが可能となる．\nさらに，本研究では，従来の線形回帰モデルでは捉えきれなかった非線形性や交互作用\n効果をスパース推定法を活用することで反映させることができる．非線形回帰や機械学習\n技術と組み合わせることにより，より精緻な価格形成モデルが構築され，従来の方法では\n見過ごされがちな要因を捉えることが可能となる．例えば，不動産市場における立地や面\n積の相互作用，消費財におけるブランド力と機能性の複合的影響など，これらの複雑な関\n係性を反映させることで，より現実的な市場分析が行えるようになる．\nまた，こうした高度な統計手法や機械学習技術の活用により，従来のヘドニック・アプ\nローチを大きく強化することができる．特に，従来のアプローチでは捉えきれなかった価\n格形成における微細な要因を，スパース推定を用いることで明示化し，より高精度な価格\n評価を可能にする．これにより，不動産市場や消費財市場における価格形成メカニズムの\n理解が深まり，実務的な応用に貢献することが期待される．\n*使用するファイル [#pe64d1b6]\n本研究で使用するファイルはすべてGoogleDriveにある。\n&br;\nuser:iie.lab.tpu.2324@gmail.com\n&br;\n\n「マイドライブ」→「学生」→「24_o3中島」→「保存期間5年」→「３）卒論（プログラム）」の中にある。\n&br;\n「引継ぎ用.zip」をダウンロードする。\n&br;\nファイルの中身は以下のようになっている。\nこれを上から順に実行していく過程を以下の章に示す。\n\n|フォルダ名|用途|ファイル名|ファイルの場所|\n|データセットの作成|犯罪データの取得|get_crime_data.py|データセットの作成\get_crime_data.py|\n|データセットの作成|統計データなどの作成|hop_統計.ipynb|データセットの作成\hop_統計.ipynb|\n|データセットの作成|取引データのメッシュ変換|緯度経度.py|データセットの作成\緯度経度.py|\n|データセットの作成|各データを取引データにまとめる|hop_メッシュ.ipynb|データセットの作成\hop_メッシュ.ipynb|\n|データセットの作成|犯罪場所との最短距離|distance.ipynb|データセットの作成\distance.ipynb|\n|前処理|データセットの前処理|前処理.ipynb|前処理\前処理.ipynb|\n|変数選択|変数の選択とハイパーパラメータのチューニング|変数選択__削除.ipynb|変数選択\変数選択_削除.ipynb|\n|要因分析|選択された変数を用いたヘドニック・アプローチによる要因分析|hedo_削除.ipynb|要因分析\hedo_削除.ipynb|\n~\nコードの変更点はディレクトリを変更するだけで大丈夫.\n\npythonのバージョンは12.0.0を用いる.\n\n使用したモジュールのバージョンは下になります。（すべては使用しておりません）\n#ref(モジュール１.png,,)\n#ref(モジュール2.png,,)\n#ref(モジュール3.png,,)\n#ref(モジュール4.png,,)\n** 犯罪発生データの取得 [#k37ae06a]\n\n「''get_crime_data.py''」を実行すると、富山県警が公開している「犯罪発生マップ」（http://www.machi-info.jp/machikado/police_pref_toyama/index.jsp）から犯罪発生データを取得し、「''crimes.csv''」として保存されます。\n&br;\n\n**建築年数の算出 [#i8ed6a27]\n建築年の年号を西暦にするために以下の式を書いた列を挿入する。\n&br;\n =IF(LEFT(M2,2)="令和", 2018 + MID(M2,3,LEN(M2)-3),\n IF(LEFT(M2,2)="平成", 1988 + MID(M2,3,LEN(M2)-3),\n IF(LEFT(M2,2)="昭和", 1925 + MID(M2,3,LEN(M2)-3), "エラー")))\n次に取引時点から西暦を抜き出す。\n&br;\n =LEFT(W2,4)\nそして、取引時点の西暦から建築年の西暦を引くことで建築年数を算出する。\n\n**使用するデータ\n実際の取引データはhttps://www.reinfolib.mlit.go.jp/realEstatePrices/ このサイトからダウンロードする。最新版を用いる場合はこちらからダウンロードを。\n\nデータセットの作成\data\estats_mesh、データセットの作成\data\estats_blockのデータはhttps://www.e-stat.go.jp/から取得できる。最新版を取得する場合はこちらから。\nhttps://www.e-stat.go.jp/gis/statmap-search?\npage=1&type=1&toukeiCode=00200521&toukeiYear=2020&aggregateUnit=H&serveyId=H002005112020&statsId=T001108&datum=2000&prefCode=16\nhttps://www.e-stat.go.jp/gis/statmap-search?page=1&type=1&toukeiCode=00200521&toukeiYear=2020&aggregateUnit=H&serveyId=H002005112020&statsId=T001101&datum=2000&prefCode=16 この二つのurlからメッシュごとのデータ。\nhttps://www.e-stat.go.jp/stat-search/files?page=1&layout=datalist&toukei=00200521&tstat=000001136464&cycle=0&tclass1=000001136472&tclass2=000001159889&cycle_facet=tclass1%3Acycle&tclass3val=0 このurlから章地域ごとのデータを。\n** データの作成と収集（統計データや施設データ）[#x2b1cdcf]\n「''hop_統計.ipynb''」を実行すると、取引データのメッシュごとの統計データの作成が行われる.\nNAVITIMEからのスクレイピングや地図画像からのデータ作成を行い、メッシュごとのデータを最終的にdata_by_mesh.csvに格納する.\n\n&br;\n df_torihiki = pd.read_csv(\n HERE / "df_housing_toyamacityメッシュ入り.csv",\n index_col = "住所",\n encoding = 'cp932') \n\nメッシュコードを文字列として格納 \n&br;\n meshes = df_torihiki["mesh_code"].astype(str).tolist()\n作成する統計データは以下のようになる。\n#ref(作成する統計データ.png,,)\n&br;\n import pandas as pd\n # ここで上で指定した説明変数を生成\n\n df_estats_block = pd.read_csv(\n DATA_DIR / "estats_block.csv",\n index_col = "KEY_CODE")\n df_estats_mesh = pd.read_csv(\n DATA_DIR / "estats_mesh.csv",\n index_col = "KEY_CODE")\n\n df_selected_block = pd.DataFrame()\n df_selected_mesh = pd.DataFrame()\n\n for line in stats_data_text.splitlines():\n \n formula, label = line.rsplit(", ", 2)\n \n if len( formula.split(" ") ) == 1:\n \n if formula in df_estats_block.columns:\n \n df_selected_block.loc[:, label] = df_estats_block.loc[:, formula]\n \n if formula in df_estats_mesh.columns:\n \n df_selected_mesh.loc[:, label] = df_estats_mesh.loc[:, formula]\n \n elif len( formula.split(" ") ) >= 2:\n \n contains_keys = re.findall("T.........", formula)\n \n if contains_keys[0] in df_estats_block.columns:\n \n for key in contains_keys:\n \n formula = formula[:formula.find(key)] + "df_estats_block.loc[:, '" + \\n formula[formula.find(key) : formula.find(key) + len(key)] + \\n "']" + formula[formula.find(key) + len(key):]\n \n exec(f"df_selected_block.loc[:, '{label}'] = {formula}")\n \n elif contains_keys[0] in df_estats_mesh.columns:\n \n for key in contains_keys:\n \n formula = formula[:formula.find(key)] + "df_estats_mesh.loc[:, '" + \\n formula[formula.find(key) : formula.find(key) + len(key)] + \\n "']" + formula[formula.find(key) + len(key):]\n \n exec(f"df_selected_mesh.loc[:, '{label}'] = {formula}")\n \n df_selected_mesh = df_selected_mesh.fillna(0)\n df_selected_block = df_selected_block.fillna(0)\n df_selected_mesh\n df_selected_block\n&br;\n # 次は施設データ．\n # 統計データと同じで，カンマの左側の数字はNAVITIMEのカテゴリID\n # コンビニ一覧のURLは "https://www.navitime.co.jp/category/0201/" なので，カテゴ \n リIDは「0201」\n\n building_data_text = '''\\n 0201, コンビニ\n 0802001, 駅\n 0805, 駐車場\n 0812, 駐輪場\n 0501, 金融機関\n 0603, 旅館\n 0608, ホテル\n 0202, スーパー\n 0205, ショッピングモール\n 0204, デパート\n 0502002, 警察署／交番\n 0504001, 小学校\n 0504002, 中学校\n 0504003, 高校\n 0504004, 大学\n 0504006, 保育園／幼稚園\n 0101, レジャー施設'''\n&br;\n # ここで，NAVITIMEからスクレイピング\n\n df_buildings = pd.DataFrame()\n\n def address2coords(address: str):\n \n for trials in range(5):\n \n res = requests.get(\n urllib.parse.quote(\n string = f"https://msearch.gsi.go.jp/address-search/AddressSearch?q={address}",\n safe = "://?="\n )\n )\n \n if res.status_code == 200:\n \n break\n \n else:\n \n time.sleep(5 * trials)\n \n else:\n \n raise Exception("国土地理院のジオコーディングAPIに接続できません。")\n \n result = res.json()\n \n time.sleep(0.5)\n # サーバに負担をかけないように0.5秒待つ\n \n if len(result) == 0:\n \n return None\n \n else:\n \n return result[0]["geometry"]["coordinates"]\n # 複数ある場合でも，0番目を採用\n\n def get_building_info(category):\n \n df_buildings = pd.DataFrame(columns = ["name", "address", "lng", "lat"])\n\n for i in range(1, sys.maxsize):\n # "while True"と同義\n \n if i == 1:\n \n url = f"https://www.navitime.co.jp/category/{category}/16201/"\n \n else:\n \n url = f"https://www.navitime.co.jp/category/{category}/16201/?page={i}"\n \n for trials in range(5):\n \n res = requests.get(url)\n \n if res.status_code == 200:\n \n break\n \n else:\n \n time.sleep(5)\n \n else:\n \n break\n # NAVITIME から 5 回連続で取得できなかったら、終了\n \n soup = BeautifulSoup(res.text, "html.parser")\n \n buildings = []\n for element in soup.find_all(class_ = "spot-name-text"):\n \n content = element.get_text()\n buildings += [{"name": content}]\n \n count = 0\n for element in soup.find_all(class_ = "spot-detail-value-text"):\n \n content = element.get_text()\n \n matches = re.findall(r"(北海道|東京都|大阪府|京都府|神奈川県|和歌山県|鹿児島県|..県)(.{1,6}[市郡区町村])", content)\n # 完璧ではないが、これで住所かを判定\n \n if len(matches) > 0:\n \n address = content\n coords = address2coords(address)\n print(address, coords)\n \n buildings[count]["address"] = address\n \n if coords == None:\n buildings[count]["lng"] = None\n buildings[count]["lat"] = None\n else:\n buildings[count]["lng"] = coords[0]\n buildings[count]["lat"] = coords[1]\n \n count += 1\n \n for building in buildings:\n \n df_buildings = add_dict_to_df(building, df_buildings)\n \n return df_buildings\n\n for line in building_data_text.splitlines():\n \n code, label = line.rsplit(", ", 2)\n \n if not (DATA_DIR / f"{code}.csv").exists():\n \n df_this_buildings = get_building_info(code)\n df_this_buildings.to_csv(DATA_DIR / f"{code}.csv", index = False)\n \n else:\n \n df_this_buildings = pd.read_csv(DATA_DIR / f"{code}.csv")\n \n df_buildings[label] = None\n \n for index, row in df_this_buildings.iterrows():\n \n if not index in df_buildings.index:\n \n df_buildings.loc[index] = pd.Series(dtype='float64')\n \n df_buildings.at[index, label] = f"{row['lng']},{row['lat']}"\n\n&br;\n df_buildings\n\n&br;\n # ここでグリッドセルごとのデータを作成する\n # さっき作った統計データ：そのまま入れる\n # さっき作った施設データ：そのグリッドセルにある数と，そのグリッドセルまでの \n 最短距離を求め，入れる\n # 地図画像にもとづく特徴量：\n # MapboxのAPIを使って，地図画像を取得する．\n # 特定のRGB値（たとえば，建物だったら[255, 0, 0]）のピクセル数を求め，画像 \n サイズの512 * 512で割って（面積比率），入れる．\n # さらに，道路だけを抜き出した2値画像をつくり，道路ネットワークを抽出．ノ \n ード数，エッジ数，密度，平均次数を求めて，入れる．\n\n df_data_by_mesh = pd.DataFrame(columns = ["mesh_code"])\n df_data_by_mesh = df_data_by_mesh.set_index("mesh_code")\n\n _count = 1\n for meshcode, coords in meshes.items():\n \n outputs = get_data_in_this_mesh(meshcode, coords)\n df_data_by_mesh = add_dict_to_df(outputs["data"], df_data_by_mesh, \n outputs["meshcode"])\n \n print(f"{_count} / {len(meshes)}")\n _count += 1\n取引データに存在するメッシュごとの統計データを作成する。\n&br;\n df_data_by_mesh\n df_data_by_mesh.to_csv(HERE / "data_by_mesh.csv")\n \n**緯度経度からメッシュを格納 [#kea42b7d]\n\n「''緯度経度.py''」を実行すると、地理情報（緯度・経度）から日本独自のメッシュコードを計算し、住宅データなどのCSVファイルにそのコードを追加し、「''df_housing_toyamacityメッシュ入り.csv''」を生成します。\n&br;\n import pandas as pd\n def latlon_to_mesh(lat, lon):\n """\n 緯度経度から9桁のメッシュコード（2分の1地域メッシュ）を取得する。\n :param lat: 緯度（度）\n :param lon: 経度（度）\n :return: 9桁のメッシュコード（文字列）\n """\n # 1次メッシュ\n mesh1_lat = int(lat * 1.5)\n mesh1_lon = int(lon - 100)\n # 2次メッシュ\n mesh2_lat = int((lat * 1.5 - mesh1_lat) * 8)\n mesh2_lon = int((lon - 100 - mesh1_lon) * 8)\n # 3次メッシュ\n mesh3_lat = int(((lat * 1.5 - mesh1_lat) * 8 - mesh2_lat) * 10)\n mesh3_lon = int(((lon - 100 - mesh1_lon) * 8 - mesh2_lon) * 10)\n # 4次メッシュ（2分の1地域メッシュ）\n mesh4_lat = int((((lat * 1.5 - mesh1_lat) * 8 - mesh2_lat) * 10 - mesh3_lat) * 2)\n mesh4_lon = int((((lon - 100 - mesh1_lon) * 8 - mesh2_lon) * 10 - mesh3_lon) * 2)\n # メッシュコードの組み立て（9桁になるよう修正）\n mesh_code = f"{mesh1_lat:01}{mesh1_lon:01}{mesh2_lat:01}{mesh2_lon:01}{mesh3_lat:01}{mesh3_lon:01}{mesh4_lat:01}{mesh4_lon:01}"\n # 最後の2桁のみ取り出して、9桁にする\n mesh_code = mesh_code[:9]\n return mesh_code\n # CSVファイルの読み込み\n file_path = r"C:\Users\nt011\Desktop\研究\富山取引のほう富山市のみ \n \df_housing_toyamacity.csv"\n df = pd.read_csv(file_path, encoding='utf-8')\n # メッシュコードの計算\n df["mesh_code"] = df.apply(lambda row: latlon_to_mesh(row["緯度"], row["経度"]), axis=1)\n cols = ["mesh_code"] + [col for col in df.columns if col != "mesh_code"]\n df = df[cols]\n # 結果のCSVを保存\n output_path = r"C:\Users\nt011\Desktop\研究\富山取引のほう富山市のみ \n \df_housing_toyamacityメッシュ入り.csv"\n df.to_csv(output_path, index=False)\n print(f"処理完了: {output_path}")\n#ref(作成されるデータ.png,,)\n\n\n\n「''hop_メッシュ.ipynb''」を実行すると、メッシュコードごとにデータを結合し「''取引_統計データ付き.csv''」に生成します。緯度・経度をもとに、500mメッシュ単位で犯罪データを分類できます。\n # メッシュGeoJsonを，分かりやすい形に直す\n _meshes = {}\n for mesh in meshes:\n meshcode = mesh["properties"]["code"]\n coords = mesh["geometry"]["coordinates"][0]\n north = coords[1][1]\n south = coords[0][1]\n east = coords[2][0]\n west = coords[0][0]\n center = [(east + west) / 2, (north + south) / 2]\n _meshes[meshcode] = {\n "north": north,\n "south": south,\n "east": east,\n "west": west,\n "center": center\n }\n meshes = _meshes\n\n # 犯罪データの座標から，どのグリッドセルに含まれるかを求め\n # 「mesh_code」列に格納する\n # 読み込んだメッシュGeoJsonに含まれない犯罪データは除外\n df_temp = df_crimes[0:0].copy()\n df_temp["mesh_code"] = None\n for index, row in df_crimes.iterrows():\n for meshcode, coords in meshes.items():\n if ((coords["south"] <= row["latitude"] < coords["north"]) and\n (coords["west"] <= row["longitude"] < coords["east"])):\n if not pd.isnull(row["datetime"]):\n df_temp.loc[index] = pd.concat([\n row,\n pd.Series({\n "mesh_code": meshcode\n })\n ])\n break\n print(f"{df_crimes.index.get_loc(index)} / {len(df_crimes)}")\n df_crimes = df_temp.copy()\n df_crimes["datetime"] = pd.to_datetime(df_crimes["datetime"], format = "%Y/%m/%d %H:%M")\n\n\n # CSVファイルの読み込み\n df_crimes = pd.read_csv(HERE / "犯罪_メッシュ入り.csv")\n #メッシュごとに件数をカウント\n mesh_counts = df_crimes['mesh_code'].value_counts().reset_index()\n mesh_counts.columns = ['mesh_code', '件数']\n # 結果を表示\n print(mesh_counts)\n # 結果をCSVファイルに保存する場合\n mesh_counts.to_csv(HERE / "犯罪_メッシュ件数.csv", index=False, encoding='utf-8-sig')\n\n&br;\n # 件数データ\n df_number = pd.read_csv(\n HERE / "犯罪_メッシュ件数.csv",\n index_col = "mesh_code")\n\n&br;\n # 取引データ\n df_torihiki = pd.read_csv(\n HERE / "df_housing_toyamacityメッシュ入り.csv",\n index_col = "mesh_code",\n encoding = 'cp932')\n\n&br;\n # メッシュコードをキーとして結合\n df_merged = pd.merge(df_torihiki, df_number, on='mesh_code', how='left', \n suffixes=('', '_number'))\n\n\n print(df_merged.head())\n # 結果をCSVファイルに保存する場合\n df_merged.to_csv(HERE / "取引_件数付き.csv", encoding='utf-8-sig')\n \n&br;\n df_mesh = pd.read_csv(\n HERE / "data_by_mesh.csv",\n index_col = "mesh_code",\n encoding = "cp932")\n\n&br;\n df_torihiki2 = pd.read_csv(\n HERE / "取引_件数付き.csv",\n index_col = "mesh_code")\n\n&br;\n # メッシュコードでデータを結合\n merged_df = pd.merge(df_torihiki2, df_mesh, on='mesh_code', how='left')\n print(merged_df.head())\n merged_df.to_csv(HERE / "取引_統計データ付き.csv", index=False, \n encoding='utf-8-sig')\n\n犯罪発生率(件/250000)はhop_メッシュ.ipynbを行った後にcsvファイル上の新しい列(犯罪発生率(件/250000)を追加して行う。\n計算方法は=件数/250000で行う。\n\n#ref(統計データ.png,,20%)\n\n「''distance.ipynb''」を実行すると、住宅地データ (取引_統計データ付き.csv) と犯罪データ (犯罪_メッシュ入り.csv) を読み込み、各住宅地から最寄りの犯罪発生地点までの最短距離を計算するために、geopy.distanceモジュールのgeodesic関数を使用します。この関数は、2点間の地球上での距離（メートル単位）を計算します。契約された土地から最短の距離で犯罪が起きた場所までの最短距離を「''取引_統計データ付き.csv''」に追加した「''住宅地_犯罪_最短距離.csv''」を生成します。\n import pandas as pd\n from geopy.distance import geodesic\n # ファイルを読み込む\n df_housing = pd.read_csv("C:/Users/nt011/Desktop/統計正しいほう/データセット \n の作成/取引_統計データ付き.csv")\n df_crimes = pd.read_csv("C:/Users/nt011/Desktop/統計正しいほう/データセット \n の作成/犯罪_メッシュ入り.csv")\n print(df_housing.columns)\n # 1. 住宅地データの緯度経度情報を抽出\n df_housing['latitude'] = df_housing['緯度']\n df_housing['longitude'] = df_housing['経度']\n # 2. 犯罪データの緯度経度情報を抽出\n df_crimes['latitude'] = df_crimes['latitude']\n df_crimes['longitude'] = df_crimes['longitude']\n # 3.最短距離を計算する関数\n def calculate_min_distance(lat1, lon1, df_target, lat_col, lon_col):\n distances = df_target.apply(lambda row: geodesic((lat1, lon1), \n (row[lat_col], row[lon_col])).meters, axis=1)\n return distances.min()\n # 4. 各住宅地から犯罪発生場所までの最短距離を計算\n df_housing['min_distance_to_crime'] = df_housing.apply(lambda row: \n calculate_min_distance(row['latitude'], row['longitude'], df_crimes, 'latitude', 'longitude'), axis=1)\n # 距離が0になる住宅地の行を表示\n df_housing[df_housing['min_distance_to_crime'] == 0]\n # 必要なカラムを選択して新しいDataFrameに保存\n df_result = df_housing[['mesh_code', '住所', '経度', '緯度', '最寄駅：名称', '最寄駅：距離（分）', '取引価格（総額）', '面積（㎡）', '取引価格（㎡単価）', \n '土地の形状', '間口', '延床面積（㎡）','建築年数', '建物の構造', \n '用途', '前面道路：方位', '前面道路：種類', '前面道路：幅員（ｍ）', '都市計画', '建ぺい率（％）', \n '容積率（％）', '取引時点', '件数', '犯罪発生率(件/250000)', '人口', '18歳未満人口割合', \n '65歳以上人口割合', '外国人人口割合', '世帯数', '単身世帯割合', '核家族世帯割合', '正規労働者割合', \n '非正規労働者割合', '最終学歴が中学以下の人口割合', '最終学歴が高校の人口割合', '最終学歴が大学以上の人口割合', \n '居住年数5年未満の人口割合', '居住年数が20年以上の人口割合', '一戸建て割合', 'アパート・低中層マンション割合', \n '高層マンション割合', 'コンビニ', '駅', '駐車場', '駐輪場', '金融機関', '旅館', 'ホテル', 'スーパー', 'ショッピングモール', 'デパート', '警察署／交番', '小学校', '中学校', '高校', '大学',\n '保育園／幼稚園', 'レジャー施設','道路の面積比率', '建物の面積比率', '水の面積比率', '空き地の面積比率', \n 'ノード数（道路）', 'エッジ数（道路）', '密度（道路）', '平均次数（道路）', 'min_distance_to_crime']]\n # 結果の保存\n df_result.to_csv('C:/Users/nt011/Desktop/統計正しいほう/データセットの作成/ \n 住宅地_犯罪_最短距離.csv', index=False, encoding='utf-8-sig')\n\n#ref(最短.png,,20%)&br;\n&br;\n** データセットの前処理 [#u50edad0]\n\n作成したデータセットを回帰分析がかけられる形に整える。\n&br;\n\n「''前処理.ipynb''」を実行し、「''住宅地_犯罪_最短距離.csv''」から「建築年数」と「最寄駅：距離（分）」の列に関する欠損値を含む行を削除した「''住宅地_犯罪_最短距離_欠損値削除後.csv''」を生成します。\n&br;\n import pandas as pd\n # CSVファイルのパス\n file_path = r"C:\Users\nt011\Desktop\統計正しいほう\前処理\住宅地_犯罪_最短距離.csv"\n # CSVファイルを読み込む\n df = pd.read_csv(file_path)\n # 欠損値を含む行を削除\n df = df.dropna()\n # 「建築年数」列に `-1` または `"VALUE"` を含む行を削除\n df = df[~df["建築年数"].astype(str).isin(["-1", "#VALUE!"])]\n # 「最寄駅：距離（分）」の値を変換\n replace_dict = {\n "30分~60分": 45,\n "1H~1H30": 75,\n "1H30~2H": 105\n }\n df["最寄駅：距離（分）"] = df["最寄駅：距離（分）"].replace(replace_dict)\n # クリーンなデータを新しいCSVファイルとして保存\n output_path = r"C:\Users\nt011\Desktop\統計正しいほう\前処理\住宅地_犯罪_最短距離_欠損値削除後.csv"\n df.to_csv(output_path, index=False, encoding='utf-8-sig')\n print("データ処理完了: 欠損値と指定された異常値を削除し、駅距離の値を変換しました。")\n print(f"新しいCSVファイルが {output_path} に保存されました。")\n\n\n生成した「''住宅地_犯罪_最短距離_欠損値削除後.csv''」から連続変数の外れ値を検出するカラムリストを定義し、データ型を数値に変換したのち検出された列を削除したリスト「''外れ値削除.csv''」を生成します。ダミー変数の処理を行う前にcsvファイル内の建物の構造の英語の部分を変換する必要がある。\n\n◆不動産取引価格~\nhttps://www.reinfolib.mlit.go.jp/realEstatePrices/~\n~\nこのurlの制度と用語のページを参考にして変換を行う。\n#ref(建物の構造.png,,30%)\n&br;\n import pandas as pd\n # CSVファイルのパス\n file_path = "C:/Users/nt011/Desktop/統計正しいほう/前処理/住宅地_犯罪_最短距 \n 離_欠損値削除後.csv"\n # CSVファイルを読み込む\n df = pd.read_csv(file_path)\n # 外れ値を検出する連続変数のカラムリスト\n continuous_columns = [\n "最寄駅：距離（分）", "取引価格（総額）", "面積（㎡）", "取引価格（㎡単価）",\n "間口", "延床面積（㎡）", "建築年数", "前面道路：幅員（ｍ）",\n "建ぺい率（％）", "容積率（％）", "犯罪発生率(件/250000)", "人口",\n "18歳未満人口割合", "65歳以上人口割合", "外国人人口割合", "世帯数",\n "単身世帯割合", "核家族世帯割合", "正規労働者割合", "非正規労働者割合",\n "最終学歴が中学以下の人口割合", "最終学歴が高校の人口割合", "最終学歴が大学以上の人口割合",\n "居住年数5年未満の人口割合", "居住年数が20年以上の人口割合", "一戸建て割合",\n "アパート・低中層マンション割合", "高層マンション割合", "道路の面積比率",\n "建物の面積比率", "水の面積比率", "空き地の面積比率", "ノード数（道路）",\n "エッジ数（道路）", "密度（道路）", "平均次数（道路）", "min_distance_to_crime"\n ]\n # データ型を数値に変換（エラーが出たらNaNにする）\n for column in continuous_columns:\n if column in df.columns:\n df[column] = pd.to_numeric(df[column], errors='coerce')\n # 外れ値を持つ行のインデックスを格納するリスト\n outlier_indices = set()\n # 外れ値の検出と削除\n for column in continuous_columns:\n if column in df.columns:\n # 欠損値を除外\n df_clean = df[column].dropna()\n # IQR（四分位範囲）の計算\n Q1 = df_clean.quantile(0.25)\n Q3 = df_clean.quantile(0.75)\n IQR = Q3 - Q1\n # 外れ値の閾値を設定\n lower_bound = Q1 - 1.5 * IQR\n upper_bound = Q3 + 1.5 * IQR\n # 外れ値のインデックスを取得\n outliers = df[(df[column] < lower_bound) | (df[column] > upper_bound)]\n outlier_indices.update(outliers.index)\n # 外れ値の行を削除\n df_cleaned = df.drop(index=outlier_indices)\n # 結果を新しいCSVに保存\n output_path = "C:/Users/nt011/Desktop/統計正しいほう/前処理/外れ値削除.csv"\n df_cleaned.to_csv(output_path, encoding="utf-8-sig", index=False)\n print(f"外れ値を削除したデータを {output_path} に保存しました。")\n print(f"削除した行数: {len(outlier_indices)}")\n print(f"最終データの行数: {df_cleaned.shape[0]}")\n\n\nそれぞれのカラムの要素を必要に応じてダミー変数化した後、データフレームを保存した「''削除_ダミー後.csv''」を生成します。~\n import pandas as pd\n # CSVファイルを読み込み\n df_data = pd.read_csv("C:/Users/nt011/Desktop/統計正しいほう/前処理/外れ値削除.csv")\n # カラム名を取得\n columns = df_data.columns.tolist()\n # 建物の構造列に「RC」がある場合、「鉄筋コンクリート造」に変換\n df_data['建物の構造'] = df_data['建物の構造'].replace('RC', '鉄筋コンクリート造')\n # 最寄り駅名を二値化（列名を「最寄り駅名_富山」に変更）\n index = columns.index("最寄駅：名称") # 変換前の位置を取得\n df_data.insert(index + 1, '最寄り駅名_富山', df_data['最寄駅：名称'].apply(lambda x: 1 if x == '富山' else 0))\n df_data = df_data.drop(columns=['最寄駅：名称']) # 元の列を削除\n # 都市計画を二値化（列名を「都市計画_市街化調整区域」に変更）\n index = columns.index("都市計画")\n df_data.insert(index + 1, '都市計画_市街化調整区域', df_data['都市計画'].apply(lambda x: 1 if x == '市街化調整区域' else 0))\n df_data = df_data.drop(columns=['都市計画'])\n # 土地の形状を二値化（列名を「土地の形状_不整形」に変更）\n index = columns.index("土地の形状")\n df_data.insert(index + 1, '土地の形状_不整形', df_data['土地の形状'].apply(lambda x: 1 if x == '不整形' else 0))\n df_data = df_data.drop(columns=['土地の形状'])\n # 用途を二値化（列名を「用途_住宅」に変更）\n index = columns.index("用途")\n df_data.insert(index + 1, '用途_住宅', df_data['用途'].apply(lambda x: 1 if x == '住宅' else 0))\n df_data = df_data.drop(columns=['用途'])\n # 建物の構造を二値化（列名を「建物の構造_木造」に変更）\n index = columns.index("建物の構造")\n df_data.insert(index + 1, '建物の構造_木造', df_data['建物の構造'].apply(lambda x: 1 if x == '木造' else 0))\n df_data = df_data.drop(columns=['建物の構造'])\n # 「前面道路：方位」のダミー変数化（基準値「北西」を0に設定）\n index = columns.index("前面道路：方位") # 変換前の位置を取得\n dummies = pd.get_dummies(df_data['前面道路：方位'], prefix='前面道路：方位', dtype=int)\n # 基準値は削除\n baseline_value = '北西'\n if f'前面道路：方位_{baseline_value}' in dummies.columns:\n dummies = dummies.drop(columns=[f'前面道路：方位_{baseline_value}'])\n # 挿入（前面道路：方位のダミー変数）\n ordered_dummies = ['前面道路：方位_西', '前面道路：方位_東', '前面道路：方位_南西', '前面道路：方位_南東', \n '前面道路：方位_南', '前面道路：方位_北東', '前面道路：方位_北']\n # 順番に従って挿入\n for col in ordered_dummies:\n if col in dummies.columns:\n df_data.insert(index + 1, col, dummies[col])\n # 「前面道路：種類」のダミー変数化\n road_types = ['市道', '私道', '県道']\n # 順番を指定してダミー変数を挿入\n for road in road_types:\n df_data.insert(index + 1, f'前面道路：種類_{road}', (df_data['前面道路：種類'] == road).astype(int))\n # 元の「前面道路：方位」および「前面道路：種類」列を削除\n df_data = df_data.drop(columns=['前面道路：方位', '前面道路：種類'])\n # 確認用にダミー変数化後のデータフレームを表示\n print(df_data.head())\n # ダミー変数化後のデータを保存\n df_data.to_csv("C:/Users/nt011/Desktop/統計正しいほう/前処理/削除_ダミー後.csv", \n index=False, encoding='utf-8-sig')\n print("ダミー変数変換完了！データが保存されました。")\n\n#ref(ダミー.png,,20%)\n&br;\n\n連続変数とダミー変数を手動で分類し、連続変数のみ二次項の作成をする。\n交互作用項の作成したデータフレームを保存した「''削除_変数作成後.csv''」を生成します。\n import pandas as pd\n import numpy as np\n from itertools import combinations\n # CSVファイルのパスを指定\n file_path = "C:/Users/nt011/Desktop/統計正しいほう/前処理/削除_ダミー後.csv"\n # UTF-8で読み込む\n df = pd.read_csv(file_path)\n # データの先頭を表示\n print(df.head())\n # 連続変数とダミー変数を手動で分類\n continuous_vars = [\n '最寄駅：距離（分）', '面積（㎡）', '間口', '延床面積（㎡）', '建築年数', '前面道路：幅員（ｍ）',\n '建ぺい率（％）', '容積率（％）', '犯罪発生率(件/250000)', '人口', '18歳未満人口割合',\n '65歳以上人口割合', '外国人人口割合', '世帯数', '単身世帯割合', '核家族世帯割合', \n '正規労働者割合', '非正規労働者割合', '最終学歴が中学以下の人口割合', '最終学歴が高校の人口割合',\n '最終学歴が大学以上の人口割合', '居住年数5年未満の人口割合', '居住年数が20年以上の人口割合',\n '一戸建て割合', 'アパート・低中層マンション割合', 'コンビニ', '駅', '駐車場', '駐輪場', \n '金融機関', '旅館', 'ホテル', 'スーパー', 'ショッピングモール', \n 'デパート', '警察署／交番', '小学校', '中学校', '高校', '大学', '保育園／幼稚園', \n 'レジャー施設','道路の面積比率', '建物の面積比率', '水の面積比率', '空き地の面積比率', \n 'ノード数（道路）', 'エッジ数（道路）', '密度（道路）', '平均次数（道路）', 'min_distance_to_crime'\n ]\n # ダミー変数をリストとして定義\n dummy_vars = [\n '最寄り駅名_富山','土地の形状_不整形', '建物の構造_木造', \n '用途_住宅', '前面道路：方位_北', '前面道路：方位_北東', \n '前面道路：方位_南', '前面道路：方位_南東', \n '前面道路：方位_南西', '前面道路：方位_東', '前面道路：方位_西', \n '前面道路：種類_市道', '前面道路：種類_県道', '前面道路：種類_私道', \n '都市計画_市街化調整区域'\n ]\n # 二次項の作成（連続変数のみ）\n for var in continuous_vars:\n df[f'{var}^2'] = df[var] ** 2\n # 交互作用項の作成\n # 連続変数 × 連続変数\n for var1, var2 in combinations(continuous_vars, 2):\n df[f'{var1}×{var2}'] = df[var1] * df[var2]\n # ダミー変数 × 連続変数\n for dummy_var in dummy_vars:\n for cont_var in continuous_vars:\n df[f'{dummy_var}×{cont_var}'] = df[dummy_var] * df[cont_var]\n # ダミー変数 × ダミー変数\n for var1, var2 in combinations(dummy_vars, 2):\n df[f'{var1}×{var2}'] = df[var1] * df[var2]\n # CSVに保存\n output_path = "C:/Users/nt011/Desktop/統計正しいほう/前処理/削除_変数作成後.csv"\n df.to_csv(output_path, index=False, encoding='utf-8-sig')\n print(f"処理が完了しました。生成されたデータは {output_path} に保存されました。")\n#ref(交互作用項.png,,20%)\n&br;\n\n*変数選択 [#u9327c42]\n***Elastic Netとは [#d2ec2870]\nElastic Netは、L1 正則化（Lasso）と L2 正則化（Ridge）の組み合わせにより、多重共線性があるデータに対して安定した回帰分析を行う手法です。\n\nLasso（L1正則化）：L1正則化は、回帰モデルの係数の一部をゼロにすることによって、特徴量選択を行います。これにより、重要な特徴量だけを残すことができますが、多重共線性が高い場合には不安定になることがあります。\n&br;\n\nRidge（L2正則化）：L2正則化は、回帰モデルの係数を小さくすることで、過学習を防ぎます。L2正則化は係数をゼロにすることはなく、むしろ係数の大きさを均等に抑える傾向がありますが、特徴量選択の機能はありません。\n&br;\n\nStandardScaler() を用いて X（説明変数）と y（目的変数）を標準化 しているため、変数のスケールによる影響を排除し、モデルの安定性を向上させています。\n\n\nAdaptive Elastic Netを用いて土地の価格に影響する重要な特徴量の係数を取得し、影響の大きい特徴量を特定します。&br;\n#ref(式.png,,50%)\n***手順 [#d2ec2870]\n「''変数選択__削除.ipynb''」を実行すると、重要な変数のみを抽出し、「''選ばれた変数.csv''」に保存します。\n import numbers\n import warnings\n import cvxpy\n import numpy as np\n from asgl import ASGL\n from sklearn.base import MultiOutputMixin\n from sklearn.base import RegressorMixin\n from sklearn.exceptions import ConvergenceWarning\n from sklearn.linear_model import ElasticNet\n from sklearn.linear_model._coordinate_descent import _alpha_grid\n from sklearn.utils import check_X_y\n from sklearn.utils.validation import check_is_fitted\n class AdaptiveElasticNet(ASGL, ElasticNet, MultiOutputMixin, RegressorMixin):\n """\n Objective function and parameters as described with modifications\n to allow alpha1, alpha2, l1_ratio1, and l1_ratio2 customization.\n """\n def __init__(\n self,\n alpha1=0.021544346900318846, # First-stage ElasticNet alpha\n alpha2=0.0009443498043343188, # Second-stage ElasticNet alpha\n *,\n l1_ratio1=0.8889, # First-stage ElasticNet L1 ratio\n l1_ratio2=0.778, # Second-stage ElasticNet L1 ratio\n gamma=0.5, # Weight adjustment exponent\n fit_intercept=True,\n precompute=False,\n max_iter=10000,\n copy_X=True,\n solver=None,\n tol=None,\n positive=False,\n positive_tol=1e-3,\n random_state=None,\n eps_coef=1e-6,\n verbose=True\n ):\n params_asgl = dict(model="lm", penalization="asgl")\n if solver is not None:\n params_asgl["solver"] = solver\n if tol is not None:\n params_asgl["tol"] = tol\n super().__init__(**params_asgl)\n self.alpha1 = alpha1\n self.alpha2 = alpha2\n self.l1_ratio1 = l1_ratio1\n self.l1_ratio2 = l1_ratio2\n self.gamma = gamma\n self.fit_intercept = fit_intercept\n self.max_iter = max_iter\n self.precompute = precompute\n self.copy_X = copy_X\n self.positive = positive\n self.positive_tol = positive_tol\n self.random_state = random_state\n self.eps_coef = eps_coef\n self.verbose = verbose\n if not self.fit_intercept:\n raise NotImplementedError\n def fit(self, X, y, check_input=True):\n if check_input:\n X_copied = self.copy_X and self.fit_intercept\n X, y = self._validate_data(\n X,\n y,\n accept_sparse="csc",\n order="F",\n dtype=[np.float64, np.float32],\n copy=X_copied,\n multi_output=True,\n y_numeric=True,\n )\n # 第一段階の ElasticNet 実行\n enet_model_1 = self.elastic_net(self.l1_ratio1, alpha=self.alpha1)\n enet_model_1.fit(X, y)\n enet_coef = enet_model_1.coef_\n # 重みの計算\n weights = 1.0 / (np.maximum(np.abs(enet_coef), self.eps_coef) ** self.gamma)\n # 第二段階の最適化\n self.coef_, self.intercept_ = self._optimize_second_stage(X, y, weights)\n # モデル属性を格納\n self.enet_coef_ = enet_coef\n self.weights_ = weights\n return self\n def predict(self, X):\n check_is_fitted(self, ["coef_", "intercept_"])\n return super(ElasticNet, self).predict(X)\n def elastic_net(self, l1_ratio, **params):\n """\n Create an ElasticNet model with the specified parameters.\n """\n elastic_net = ElasticNet(l1_ratio=l1_ratio)\n for k, v in self.get_params().items():\n try:\n elastic_net = elastic_net.set_params(**{k: v})\n except ValueError:\n pass # Ignore parameters not supported by ElasticNet\n elastic_net = elastic_net.set_params(**params)\n return elastic_net\n def _optimize_second_stage(self, X, y, weights):\n """\n Perform second-stage optimization with adaptive weights.\n Returns\n -------\n coef : np.array, shape (n_features,)\n intercept : float\n """\n n_samples, n_features = X.shape\n beta_variables = [cvxpy.Variable(n_features)]\n model_prediction = 0.0\n if self.fit_intercept:\n beta_variables = [cvxpy.Variable(1)] + beta_variables\n ones = cvxpy.Constant(np.ones((n_samples, 1)))\n model_prediction += ones @ beta_variables[0]\n # モデル予測\n model_prediction += X @ beta_variables[1]\n error = cvxpy.sum_squares(y - model_prediction) / (2 * n_samples)\n l1_coefs = self.alpha2 * self.l1_ratio2\n # 第二段階の正則化項\n l1_penalty = cvxpy.Constant(l1_coefs * weights) @ cvxpy.abs(\n beta_variables[1]\n )\n l2_penalty = (\n cvxpy.Constant(self.alpha1 * (1 - self.l1_ratio1))\n * cvxpy.sum_squares(beta_variables[1])\n )\n constraints = [b >= 0 for b in beta_variables] if self.positive else []\n # 最適化問題の定義\n problem = cvxpy.Problem(\n cvxpy.Minimize(error + l1_penalty + l2_penalty), \n constraints=constraints\n )\n problem.solve(solver="OSQP", max_iter=self.max_iter)\n if problem.status != "optimal":\n raise ConvergenceWarning(\n f"Solver did not reach optimum (Status: {problem.status})"\n )\n beta_sol = np.concatenate([b.value for b in beta_variables], axis=0)\n beta_sol[np.abs(beta_sol) < self.tol] = 0\n intercept, coef = beta_sol[0], beta_sol[1:]\n coef = np.maximum(coef, 0) if self.positive else coef\n return coef, intercept\n from sklearn.linear_model import ElasticNetCV\n from sklearn.model_selection import GridSearchCV\n from sklearn.datasets import make_regression\n import numpy as np\n from sklearn.model_selection import train_test_split\n import pandas as pd\n from sklearn.preprocessing import StandardScaler\n # CSVファイルのパスを指定してください\n file_path = "C:/Users/nt011/Desktop/統計正しいほう/変数選択/削除_変数作成後.csv"\n # CSVファイルの読み込み\n df = pd.read_csv(file_path) # 文字コードが異なる場合は、'utf-8' を他のエンコーディングに変更してください\n # 特徴量とターゲットの分割\n X = df.drop(columns=['取引価格（㎡単価）']) # 取引価格（㎡単価）をyとして分離\n y = df['取引価格（㎡単価）'] # ターゲット変数\n # 🔥 標準化の実施\n scaler_X = StandardScaler()\n scaler_y = StandardScaler()\n X_scaled = scaler_X.fit_transform(X) # Xの標準化\n y_scaled = scaler_y.fit_transform(y.values.reshape(-1, 1)).ravel() # yの標準化\n # データ分割\n X_train, X_test, y_train, y_test = train_test_split(X_scaled, y_scaled, \n test_size=0.2, random_state=42)\n # 訓練データの一部を使用\n X_sample, _, y_sample, _ = train_test_split(X_train, y_train, test_size=0.8, \n random_state=42)\n # 第一段階: ElasticNetCVを使用した最適パラメータの導出\n enet_cv = ElasticNetCV(\n l1_ratio=np.linspace(0.0001, 1, 25), # l1_ratioの候補\n alphas=np.logspace(-5, 0, 25), # alphaの候補\n cv=5, # 交差検証の分割数\n random_state=42,\n n_jobs=-1,\n verbose=1\n )\n enet_cv.fit(X_train, y_train)\n # 第一段階の最適パラメータと係数を取得\n alpha1_opt = enet_cv.alpha_\n l1_ratio1_opt = enet_cv.l1_ratio_\n enet_coef = enet_cv.coef_\n print(f"第一段階の最適パラメータ: alpha1={alpha1_opt}, l1_ratio1= \n {l1_ratio1_opt}")\n\n # ベストパラメータをalpha2とl1_ratio2でそれぞれ格納\n best_alpha2 = best_params['alpha2']\n best_l1_ratio2 = best_params['l1_ratio2']\n print(best_alpha2)\n print(best_l1_ratio2)\nASGLとElasticNetを組み合わせたAdaptiveElasticNetクラスを定義し、CSVから読み込んだデータを標準化後にモデルをフィッティングして、非ゼロの係数を持つ重要な特徴量を抽出し、その結果をCSVに保存する一連の処理を実行します。\n import numbers\n import warnings\n import cvxpy\n import numpy as np\n from asgl import ASGL\n from sklearn.base import MultiOutputMixin\n from sklearn.base import RegressorMixin\n from sklearn.exceptions import ConvergenceWarning\n from sklearn.linear_model import ElasticNet\n from sklearn.linear_model._coordinate_descent import _alpha_grid\n from sklearn.utils import check_X_y\n from sklearn.utils.validation import check_is_fitted\n class AdaptiveElasticNet(ASGL, ElasticNet, MultiOutputMixin, RegressorMixin):\n """\n Objective function and parameters as described with modifications\n to allow alpha1, alpha2, l1_ratio1, and l1_ratio2 customization.\n """\n def __init__(\n self,\n alpha1=0.01333521432163324, # 第一段階で得られたalpha1を指定\n alpha2=0.034807005884284134, # Second-stage ElasticNet alpha\n *,\n l1_ratio1=1, # First-stage ElasticNet L1 ratio\n l1_ratio2=0.1667499999999998, # Second-stage ElasticNet L1 ratio\n gamma=0.5, # Weight adjustment exponent\n fit_intercept=True,\n precompute=False,\n max_iter=10000,\n copy_X=True,\n solver=None,\n tol=None,\n positive=False,\n positive_tol=1e-3,\n random_state=None,\n eps_coef=1e-6,\n verbose=True\n ):\n params_asgl = dict(model="lm", penalization="asgl")\n if solver is not None:\n params_asgl["solver"] = solver\n if tol is not None:\n params_asgl["tol"] = tol\n super().__init__(**params_asgl)\n self.alpha1 = alpha1\n self.alpha2 = alpha2\n self.l1_ratio1 = l1_ratio1\n self.l1_ratio2 = l1_ratio2\n self.gamma = gamma\n self.fit_intercept = fit_intercept\n self.max_iter = max_iter\n self.precompute = precompute\n self.copy_X = copy_X\n self.positive = positive\n self.positive_tol = positive_tol\n self.random_state = random_state\n self.eps_coef = eps_coef\n self.verbose = verbose\n if not self.fit_intercept:\n raise NotImplementedError\n def fit(self, X, y, check_input=True):\n if check_input:\n X_copied = self.copy_X and self.fit_intercept\n X, y = self._validate_data(\n X,\n y,\n accept_sparse="csc",\n order="F",\n dtype=[np.float64, np.float32],\n copy=X_copied,\n multi_output=True,\n y_numeric=True,\n )\n # 第一段階の ElasticNet 実行\n enet_model_1 = self.elastic_net(self.l1_ratio1, alpha=self.alpha1)\n enet_model_1.fit(X, y)\n enet_coef = enet_model_1.coef_\n # 重みの計算\n weights = 1.0 / (np.maximum(np.abs(enet_coef), self.eps_coef) ** self.gamma)\n # 第二段階の最適化\n self.coef_, self.intercept_ = self._optimize_second_stage(X, y, weights)\n # モデル属性を格納\n self.enet_coef_ = enet_coef\n self.weights_ = weights\n return self\n def predict(self, X):\n check_is_fitted(self, ["coef_", "intercept_"])\n return super(ElasticNet, self).predict(X)\n def elastic_net(self, l1_ratio, **params):\n """\n Create an ElasticNet model with the specified parameters.\n """\n elastic_net = ElasticNet(l1_ratio=l1_ratio)\n for k, v in self.get_params().items():\n try:\n elastic_net = elastic_net.set_params(**{k: v})\n except ValueError:\n pass # Ignore parameters not supported by ElasticNet\n elastic_net = elastic_net.set_params(**params)\n return elastic_net\n def _optimize_second_stage(self, X, y, weights):\n """\n Perform second-stage optimization with adaptive weights.\n Returns\n -------\n coef : np.array, shape (n_features,)\n intercept : float\n """\n n_samples, n_features = X.shape\n beta_variables = [cvxpy.Variable(n_features)]\n model_prediction = 0.0\n if self.fit_intercept:\n beta_variables = [cvxpy.Variable(1)] + beta_variables\n ones = cvxpy.Constant(np.ones((n_samples, 1)))\n model_prediction += ones @ beta_variables[0]\n # モデル予測\n model_prediction += X @ beta_variables[1]\n error = cvxpy.sum_squares(y - model_prediction) / (2 * n_samples)\n l1_coefs = self.alpha2 * self.l1_ratio2\n # 第二段階の正則化項\n l1_penalty = cvxpy.Constant(l1_coefs * weights) @ cvxpy.abs(\n beta_variables[1]\n )\n l2_penalty = (\n cvxpy.Constant(self.alpha1 * (1 - self.l1_ratio1))\n * cvxpy.sum_squares(beta_variables[1])\n )\n constraints = [b >= 0 for b in beta_variables] if self.positive else []\n # 最適化問題の定義\n problem = cvxpy.Problem(\n cvxpy.Minimize(error + l1_penalty + l2_penalty), \n constraints=constraints\n )\n problem.solve(solver="OSQP", max_iter=self.max_iter)\n if problem.status != "optimal":\n raise ConvergenceWarning(\n f"Solver did not reach optimum (Status: {problem.status})"\n )\n beta_sol = np.concatenate([b.value for b in beta_variables], axis=0)\n beta_sol[np.abs(beta_sol) < self.tol] = 0\n intercept, coef = beta_sol[0], beta_sol[1:]\n coef = np.maximum(coef, 0) if self.positive else coef\n return coef, intercept\n import pandas as pd\n from sklearn.preprocessing import StandardScaler\n # CSVファイルのパスを指定してください\n file_path = "C:/Users/monda/OneDrive/デスクトップ/統計正しいほう/変数選択/削除_変数作成後.csv"\n # CSVファイルの読み込み\n df = pd.read_csv(file_path) # 文字コードが異なる場合は、'utf-8' を他のエンコー ディングに変更してください\n # 特徴量とターゲットの分割\n X = df.drop(columns=['取引価格（㎡単価）']) # 取引価格（㎡単価）をyとして分離\n y = df['取引価格（㎡単価）'] # ターゲット変数\n # 🔥 標準化の実施\n scaler_X = StandardScaler()\n scaler_y = StandardScaler()\n X_scaled = scaler_X.fit_transform(X) # Xの標準化\n y_scaled = scaler_y.fit_transform(y.values.reshape(-1, 1)).ravel() # yの標準化\n model = AdaptiveElasticNet()\n model.fit(X_scaled, y_scaled)\n # 係数の表示（重要な特徴量を確認）\n coef = model.coef_\n print("特徴量の係数:", coef)\n # 係数がゼロでない特徴量を選定\n selected_features = X.columns[coef != 0]\n print("選定された特徴量:", selected_features)\n # 係数と特徴量名をデータフレームに変換\n selected_features_df = pd.DataFrame({\n 'Feature': X.columns,\n 'Coefficient': coef\n })\n # 係数がゼロでない特徴量のみをフィルタリング\n selected_features_df = selected_features_df[selected_features_df['Coefficient'] != 0]\n # CSVファイルとして保存\n output_file_path = "選ばれた変数.csv" # 保存先を指定してください\n selected_features_df.to_csv(output_file_path, index=False, encoding='utf-8-sig') \n # 日本語対応のエンコーディング\n print(f"選定された特徴量と係数を {output_file_path} に保存しました。")\n #ref(変数選択.png,,50%)\n&br;\nAENの数式は以下のようなものを扱う。そして、今回行ったハイパーパラメータは以下のものである。\nAENの推計は，二段階で行われる．まず，予備推計としてEN で係数を推定する．その\nうえで，係数の絶対値が小さい変数に対してより強い罰則を与えるようL1 ノルムの正則化\n項を変数ごとに調整したうえで，改めてEN を実施する．\n#ref(AEN.png,,150%)\n#ref(ハイパーパラメータ.png,,150%)\n\n選択された変数は以下のようになる\n#ref(選ばれた変数.png,,)\n\n*要因分析 [#w249dca6]\n***手順 [#d20580fd]\n「''hedo.ipynb''」を実行すると、有意水準に基づいて有意性を判定する関数\n import pandas as pd\n # Excelファイルを読み込み\n file_path = "c:/Users/nt011/Desktop/統計正しいほう/要因分析/位置情報_削除.csv"\n df = pd.read_csv(file_path)\n # 住所が含まれている列を特定して、住所ごとの出現回数を確認\n address_column = '住所' # 住所列の名前がわからない場合はここを変更してください\n if address_column in df.columns:\n address_counts = df[address_column].dropna().value_counts() # 住所ごとの出現回数を集計\n # 出現回数をデータフレームに変換\n address_counts_df = address_counts.reset_index()\n address_counts_df.columns = ['住所', '出現回数'] # 列名を設定\n # 結果をCSVファイルに書き出し\n output_file = "c:/Users/nt011/Desktop/統計正しいほう/要因分析/住所出現回数_削除.csv"\n address_counts_df.to_csv(output_file, index=False, encoding='utf-8-sig') \n # UTF-8エンコーディングで保存\n print(f"住所ごとの出現回数がCSVファイルに保存されました: {output_file}")\n else:\n print("住所列が見つかりません")\n回帰係数、t値、p値、有意性を格納したリストを「''回帰分析結果_削除.csv''」に生成します。\n import numpy as np\n import pandas as pd\n from sklearn.preprocessing import StandardScaler\n # ファイルパスの指定\n selected_features_path = "選ばれた変数_削除.csv"\n all_data_path = "c:/Users/nt011/Desktop/統計正しいほう/要因分析/削除_変数作成後.csv"\n location_data_path = "c:/Users/nt011/Desktop/統計正しいほう/要因分析/位置情報_削除.csv"\n output_data_path = "選択された特徴量とデータ_削除.csv"\n # 選定された変数のCSVファイルを読み込み\n selected_features_df = pd.read_csv(selected_features_path, encoding='utf-8-sig')\n selected_features = selected_features_df['Feature'].tolist()\n # 元データを再度読み込み\n all_data_df = pd.read_csv(all_data_path)\n # 選定された説明変数と目的変数を抽出\n df = all_data_df[selected_features + ['取引価格（㎡単価）']]\n # 位置情報データを読み込み\n location_df = pd.read_csv(location_data_path, usecols=['住所'])\n # 元データと位置情報を結合（共通のキーがある場合）\n df = pd.concat([df, location_df], axis=1)\n # 統合データをCSVファイルとして保存\n df.to_csv(output_data_path, index=False, encoding='utf-8-sig')\n print(f"選定された変数と位置情報を含むデータを {output_data_path} に保存しました。")\n\n import pandas as pd\n # ファイルパス\n count_file = "c:/Users/nt011/Desktop/統計正しいほう/要因分析/住所出現回数_削除.csv"\n data_file = "C:/Users/nt011/Desktop/統計正しいほう/要因分析/選択された特徴量とデータ_削除.csv"\n output_file = "C:/Users/nt011/Desktop/統計正しいほう/要因分析/ソート済みデータ_削除.csv"\n # CSVファイルの読み込み\n count_df = pd.read_csv(count_file, encoding='utf-8-sig')\n data_df = pd.read_csv(data_file, encoding='utf-8-sig')\n # 住所列の名前\n address_column = '住所'\n # 緯度・経度を除外（該当カラムがある場合）\n drop_columns = ['緯度', '経度']\n data_df = data_df.drop(columns=[col for col in drop_columns if col in \n data_df.columns])\n # 住所の出現回数をデータに統合\n merged_df = data_df.merge(count_df, on=address_column, how='left')\n # 出現回数の多い順にソート\n sorted_df = merged_df.sort_values(by='出現回数', ascending=False).drop(columns=['出現回数'])\n # ソート済みデータをCSVに保存\n sorted_df.to_csv(output_file, encoding='utf-8-sig', index=False)\n print(f"ソート済みデータを {output_file} に保存しました。")\n import pandas as pd\n import numpy as np\n from sklearn.linear_model import LinearRegression\n from sklearn.preprocessing import StandardScaler\n import statsmodels.api as sm\n # データの読み込み\n file_path = "C:/Users/nt011/Desktop/統計正しいほう/要因分析/ソート済みデータ_削除.csv"\n data = pd.read_csv(file_path)\n # 目的変数と説明変数の設定\n y = data['取引価格（㎡単価）'] # 取引価格を目的変数に\n X = data.drop(columns=['取引価格（㎡単価）', '住所']) # 取引価格と住所以外を説明変数に\n # 🔥 標準化の実施\n scaler_X = StandardScaler()\n scaler_y = StandardScaler()\n X_scaled = scaler_X.fit_transform(X) # Xの標準化\n y_scaled = scaler_y.fit_transform(y.values.reshape(-1, 1)).ravel() # yの標準化\n # 線形回帰モデルの作成\n model = LinearRegression()\n model.fit(X_scaled, y)\n # 回帰分析の結果を表示\n X_scaled_with_const = sm.add_constant(X_scaled) # 定数項を追加\n ols_model = sm.OLS(y, X_scaled_with_const).fit()\n # 必要な情報を抽出\n r_squared = ols_model.rsquared # 決定係数\n coefficients = ols_model.params # 回帰係数（定数項も含む）\n t_values = ols_model.tvalues # t値（定数項も含む）\n p_values = ols_model.pvalues # p値（有意水準を判定するため）\n variable_names = ['定数項'] + list(X.columns)\n # 有意水準に基づいて有意性を判定する関数\n def significance(p_value):\n if p_value < 0.01:\n return "1% 有意"\n elif p_value < 0.05:\n return "5% 有意"\n elif p_value < 0.10:\n return "10% 有意"\n else:\n return "有意でない"\n # 結果を格納するリスト\n results = []\n # 回帰係数、t値、p値、有意性をリストに格納\n for var_name, coef, t_val, p_val in zip(variable_names, coefficients, t_values, p_values):\n significance_level = significance(p_val)\n results.append([var_name, coef, t_val, p_val, significance_level])\n # 結果をDataFrameに変換\n results_df = pd.DataFrame(results, columns=["変数", "回帰係数", "t値", "p値", "有意性"])\n # 予測値の計算\n data["標準化予測取引価格（㎡単価）"] = ols_model.predict(X_scaled_with_const)\n data["標準化実測取引価格（㎡単価）"] = y_scaled\n # **結果を表示**\n print("\n==== 回帰分析結果 ====")\n print(f"決定係数 (R^2): {r_squared:.4f}\n")\n print(results_df)\n # 予測結果を含めたデータを保存\n predicted_output_file_path = "C:/Users/nt011/Desktop/統計正しいほう/要因分析/予測結果_削除.csv"\n data.to_csv(predicted_output_file_path, index=False, encoding="utf-8-sig")\n # CSVファイルとして保存\n output_file_path = "C:/Users/nt011/Desktop/統計正しいほう/要因分析/回帰分析結果_削除.csv"\n results_df.to_csv(output_file_path, index=False, encoding="utf-8-sig")\n print(f"回帰分析結果が{output_file_path}に保存されました。")\n#ref(ヘドニック.png,,50%)\n)
テキスト整形のルールを表示する
![[PukiWiki]](image/pukiwiki.png "[PukiWiki]")
