極座標

極座標は、平面上の一点を中心に、その点からの距離と角度で位置を表現する座標系のことを指します。

例えば、ある場所から見た時に、その場所がどれだけ離れていて、どの方向にあるかを考えてみましょう。
その場所の「距離」は、ある場所からその場所までの「直線」の長さを指し、「角度」はその場所がどの方向にあるかを表します。

この「距離」と「角度」を使って、ある場所の位置を表現するのが極座標です。
例えば、ある場所があなたの自宅から10メートル離れており、その場所があなたの自宅の右側にあるとします。
この場合、その場所の極座標は「(10メートル, 右)」と表現できます。

このように、極座標は「距離」と「角度」を使って、場所の位置を表現する座標系です。
これは、地図やコンパスなど、位置を表現するための座標系としてよく使われます。

サンプルコード

Pythonとmatplotlibを使用して複雑なグラフを作成する例として、4つの象限を持つ極座標のグラフを考えてみましょう。

この例では、Pythonのnumpyライブラリを使用してランダムなデータを生成し、そのデータを使用してグラフを描きます。

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

1
2

import numpy as np
import matplotlib.pyplot as plt

次に、ランダムなデータを生成します。

1
2

theta = np.linspace(0, 2*np.pi, 100)
r = np.random.uniform(0, 1, 100)

次に、4つの象限を持つ極座標のグラフを作成します。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19

fig, axes = plt.subplots(2, 2, subplot_kw=dict(polar=True), figsize=(16, 8))

# 第1象限
axes[0, 0].plot(theta, r)

# 第2象限
axes[0, 1].plot(theta, r)
axes[0, 1].set_theta_offset(np.pi)

# 第3象限
axes[1, 0].plot(theta, r)
axes[1, 0].set_theta_direction(-1)

# 第4象限
axes[1, 1].plot(theta, r)
axes[1, 1].set_theta_offset(np.pi)
axes[1, 1].set_theta_direction(-1)

plt.show()

このコードは、4つの象限を持つ極座標のグラフを描きます。

subplot_kw=dict(polar=True)を指定することで、各サブプロットが極座標を持つグラフになります。

set_theta_offsetメソッドとset_theta_directionメソッドを使用して、各象限の位置と方向を設定します。

なお、このコードは適当なデータを使用していますので、実際のデータを使用する場合には、そのデータを適切に準備する必要があります。

ソースコード解説

コードの概要

このコードは、numpyとmatplotlibを使用して、4つの極座標グラフを作成しています。
それぞれのグラフは極座標上の異なる象限を示しており、乱数を使って円周上に点を配置しています。

ライブラリのインポート

1
2

import numpy as np
import matplotlib.pyplot as plt

• numpyは数値計算を行うためのライブラリです。
• matplotlib.pyplotはグラフを描画するためのライブラリです。

データの準備

1
2

theta = np.linspace(0, 2*np.pi, 100)
r = np.random.uniform(0, 1, 100)

• thetaは0から2πまでの値を等間隔で100点取得しています。
• rは0から1までの乱数を100点取得しています。

グラフの作成

1

fig, axes = plt.subplots(2, 2, subplot_kw=dict(polar=True), figsize=(16, 8))

• plt.subplots()で2x2の4つの極座標グラフを作成し、figに図全体、axesにそれぞれのサブプロットが割り当てられます。

第1象限のグラフ

1

axes[0, 0].plot(theta, r)

• axes[0, 0]はsubplotの位置を指定しています。
  第1象限のグラフを作成しています。

第2象限のグラフ

1
2

axes[0, 1].plot(theta, r)
axes[0, 1].set_theta_offset(np.pi)

• axes[0, 1]は第2象限の位置を指定しています。
  set_theta_offset(np.pi)でグラフの方向を反転しています。

第3象限のグラフ

1
2

axes[1, 0].plot(theta, r)
axes[1, 0].set_theta_direction(-1)

• axes[1, 0]は第3象限の位置を指定しています。
  set_theta_direction(-1)でグラフの方向を反転しています。

第4象限のグラフ

1
2
3

axes[1, 1].plot(theta, r)
axes[1, 1].set_theta_offset(np.pi)
axes[1, 1].set_theta_direction(-1)

• axes[1, 1]は第4象限の位置を指定しています。
  set_theta_offset(np.pi)とset_theta_direction(-1)でグラフの方向を反転しています。

グラフの表示

1

plt.show()

• 作成したグラフを表示しています。
  それぞれのサブプロットが極座標上の異なる象限を示しています。

3Dプロット

Plotlyライブラリを使って、複数の3Dプロットを組み合わせて表現してみます。

以下に、複数の3Dサブプロットを表示する例を示します。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25

import plotly.graph_objs as go
import plotly.subplots as subplots
import numpy as np

# サンプルデータの生成
x = np.linspace(-5, 5, 100)
y = np.linspace(-5, 5, 100)
x, y = np.meshgrid(x, y)
z1 = np.sin(np.sqrt(x**2 + y**2))
z2 = np.cos(np.sqrt(x**2 + y**2))
z3 = x**2 - y**2

# 3Dサブプロットの作成
fig = subplots.make_subplots(rows=1, cols=3,
                             subplot_titles=('sin(sqrt(x^2 + y^2))', 'cos(sqrt(x^2 + y^2))', 'x^2 - y^2'),
                             specs=[[{'type': 'surface'}, {'type': 'surface'}, {'type': 'surface'}]])

# 3つの3Dサブプロットを追加
fig.add_trace(go.Surface(x=x, y=y, z=z1, colorscale='Viridis'), row=1, col=1)
fig.add_trace(go.Surface(x=x, y=y, z=z2, colorscale='Rainbow'), row=1, col=2)
fig.add_trace(go.Surface(x=x, y=y, z=z3, colorscale='Jet'), row=1, col=3)

# レイアウトの設定
fig.update_layout(title_text='複数の3Dサブプロット', height=600, width=1000)
fig.show()

このコードは、3つの異なる3Dプロットを1つのグラフに表示しています。

それぞれのプロットは異なる関数を表しており、それを3Dサブプロットとして横に並べて表示しています。

サンプルデータとして、2つの円周の$sin$、$cos$関数、そして$ x^2 - y^2 $の関数を生成し、それぞれの3Dサブプロットを作成しています。

これにより、複数の3Dグラフを組み合わせて表示する方法がわかります。

ソースコード解説

ソースコードの各部分を順に説明します。

ライブラリのインポート

1
2
3

import plotly.graph_objs as go
import plotly.subplots as subplots
import numpy as np

• plotly.graph_objs はPlotlyのグラフオブジェクトを作成するためのモジュールです。
• plotly.subplots はPlotlyのサブプロットを作成するためのモジュールです。
• numpy は数値計算を行うためのライブラリです。

サンプルデータの生成

1
2
3
4
5
6

x = np.linspace(-5, 5, 100)
y = np.linspace(-5, 5, 100)
x, y = np.meshgrid(x, y)
z1 = np.sin(np.sqrt(x**2 + y**2))
z2 = np.cos(np.sqrt(x**2 + y**2))
z3 = x**2 - y**2

• np.linspace() は一様に区間を分割する関数で、-5から5までの区間を100分割してxとyの値を生成しています。
• np.meshgrid() はx、yの値から格子状の座標を作成しています。
• np.sin()、np.cos()、**（べき乗演算子）を使用してz1、z2、z3の値を計算しています。

3Dサブプロットの作成

1
2
3

fig = subplots.make_subplots(rows=1, cols=3,
                             subplot_titles=('sin(sqrt(x^2 + y^2))', 'cos(sqrt(x^2 + y^2))', 'x^2 - y^2'),
                             specs=[[{'type': 'surface'}, {'type': 'surface'}, {'type': 'surface'}]])

• subplots.make_subplots() は、指定された行数と列数のサブプロットを作成する関数です。
  ここでは1行3列のサブプロットを作成しています。
• subplot_titles は、各サブプロットのタイトルを指定しています。
• specs は各サブプロットのタイプを指定しています。

3つの3Dサブプロットを追加

1
2
3

fig.add_trace(go.Surface(x=x, y=y, z=z1, colorscale='Viridis'), row=1, col=1)
fig.add_trace(go.Surface(x=x, y=y, z=z2, colorscale='Rainbow'), row=1, col=2)
fig.add_trace(go.Surface(x=x, y=y, z=z3, colorscale='Jet'), row=1, col=3)

• fig.add_trace() を使用して、各3Dサブプロットを作成しています。
• go.Surface() は3Dサーフェスプロットを作成するPlotlyの関数です。
• colorscale はカラーマップの設定を行っています。

レイアウトの設定と表示

1
2

fig.update_layout(title_text='複数の3Dサブプロット', height=600, width=1000)
fig.show()

• fig.update_layout() でグラフ全体のレイアウトを設定しています。
  ここではグラフのタイトルとサイズを指定しています。
• fig.show() で作成したグラフを表示しています。

このコードは、3つの異なる関数を表す3Dサブプロットを作成し、それらを1つのグラフにまとめて表示するものです。

それぞれのサブプロットは、異なる特徴を持つ関数を3Dプロットとして表現しています。

非線形の最適化問題

非線形の最適化問題を解くためには、SciPyライブラリのminimize()関数を使用します。

ここでは、簡単な非線形の最適化問題を解いてみます。

例として、次の非線形関数を最小化する問題を考えます。

$$
f(x) = x^2 + 5 \sin(x)
$$

この関数を最小化する$ (x) $の値を見つけることを目指します。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29

import numpy as np
from scipy.optimize import minimize
import matplotlib.pyplot as plt

# 最小化する非線形関数
def objective_function(x):
    return x**2 + 5 * np.sin(x)

# 初期値
x0 = 0

# 最適化の実行
result = minimize(objective_function, x0, method='BFGS')

# 結果の表示
print("Optimal solution:", result.x)

# グラフ化
x_vals = np.linspace(-10, 10, 100)
y_vals = objective_function(x_vals)

plt.plot(x_vals, y_vals, label='Objective Function')
plt.scatter(result.x, objective_function(result.x), color='red', label='Optimal Solution')
plt.xlabel('x')
plt.ylabel('f(x)')
plt.title('Nonlinear Optimization')
plt.legend()
plt.grid(True)
plt.show()

このコードは、minimize()関数を使用して非線形関数を最小化しています。

関数の最小値を求めるためにBFGS（Broyden-Fletcher-Goldfarb-Shanno）最適化手法を使用しています。

最適解は result.x に格納されます。

また、matplotlibライブラリを使用して、非線形関数をグラフ化しています。

最小値の位置を赤色の点で表示しています。

ソースコード解説

このコードは、PythonのSciPyライブラリを使用して非線形最適化問題を解くものです。

1. numpyおよびscipy.optimizeから必要なライブラリをインポートしています。
   また、グラフを描画するためにmatplotlib.pyplotもインポートしています。

2. objective_function(x)という関数を定義しています。
   これは最小化する非線形関数です。
   ここでは$ ( f(x) = x^2 + 5 \sin(x) ) $となっています。

3. 初期値 x0 を 0 としています。
   最適化アルゴリズムはこの初期値から始まり、関数の最小値を見つけようとします。

4. minimize()関数を使用して、objective_functionを最小化します。
   ここでは** BFGS（Broyden-Fletcher-Goldfarb-Shanno）法 **を使用しています。

5. result.xには、最適化アルゴリズムによって見つけられた最適解が格納されます。

6. グラフ化のために、np.linspace()を使用して$ ( x ) $の値を範囲$ [-10, 10] $で生成し、それに対応する$ ( f(x) ) $の値を計算しています。

7. plt.plot()を使って非線形関数を青色の折れ線グラフでプロットし、plt.scatter()を使って最適解の位置を赤い点で表示しています。

8. 最後に、グラフの軸ラベル、タイトル、凡例、グリッドを設定し、plt.show()でグラフを表示しています。

このコードは、非線形最適化問題を解く手法とその結果を可視化する方法を示しています。

関数の最小値を見つけ、最適解の位置を視覚的に示すためにグラフを使用しています。

結果解説

この非線形関数$ ( f(x) = x^2 + 5 \sin(x) ) $を最小化するために、Scipyのminimize()関数を使用しました。

最適化手法としてはBFGS法を選択しました。

結果として得られた最適解は$ ( x \approx -1.11051052 ) $です。

この$ ( x ) $の値が関数$ ( f(x) ) $を最小化する点です。

グラフでは、横軸が$ ( x )$、縦軸が$ ( f(x) ) $となっており、非線形関数が青線で表示されています。

赤い点が求められた最適解の位置を示しています。

この点は、関数が最小値を取る位置を表しています。

関数$ ( f(x) ) は ( x \approx -1.11051052 ) $の位置で最小値を持ちます。

この値は、関数が最小となる$ ( x ) $の位置を示しており、BFGS最適化アルゴリズムがこの値を見つけることができたことを示しています。

資源割り当て問題 (Resource Allocation Problem)

最適化問題として、**資源割り当て問題 (Resource Allocation Problem)**を考えてみましょう。

この問題では、限られた資源を異なるタスクに効果的に割り当て、特定の目的関数を最大化または最小化することが求められます。

例えば、以下のような状況を考えます：

問題：資源割り当て問題

資源：

ある会社のプロジェクトマネージャが5人の従業員と10,000ドルの予算を持っています。

タスク：

プロジェクトには3つの重要なタスクがあります。
それぞれのタスクには異なる従業員のスキルが必要で、予算も異なります。

目標：

全体のプロジェクト効率を最大化するために、どの従業員にどのタスクを割り当て、予算をどのように使うべきかを決定します。

この問題を数理最適化問題としてモデル化し、例えば全体のプロジェクト効率を最大化するような目的関数を設定し、従業員と予算の制約条件を考慮して解を求めることができます。

このような問題は実際のビジネスやプロジェクト管理でよく発生し、最適な資源割り当てによって企業の生産性や利益を向上させることが期待されます。

サンプルソース

資源割り当て問題を解くために、PuLPと呼ばれる線形プログラミングライブラリを使用します。

以下は、Pythonコードの一例です。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47

from pulp import LpProblem, LpVariable, lpSum, LpMaximize

def solve_resource_allocation_problem():
    # 問題の設定
    problem = LpProblem("Resource_Allocation", LpMaximize)

    # 従業員とタスクの数
    num_employees = 5
    num_tasks = 3

    # 従業員のスキル
    skills = [
        [3, 1, 4],  # 従業員1のスキル
        [2, 0, 5],  # 従業員2のスキル
        [1, 4, 2],  # 従業員3のスキル
        [4, 3, 1],  # 従業員4のスキル
        [0, 2, 3]   # 従業員5のスキル
    ]

    # タスクごとの予算
    budgets = [3000, 5000, 2000]

    # 従業員ごとの変数（0または1のバイナリ変数）
    employees = [LpVariable(f"Employee_{i}", cat='Binary') for i in range(num_employees)]

    # 目的関数（最大化する対象：全体のプロジェクト効率）
    problem += lpSum([employees[i] * skills[i][j] for i in range(num_employees) for j in range(num_tasks)]), "Total_Efficiency"

    # 制約条件
    for j in range(num_tasks):
        problem += lpSum([employees[i] * skills[i][j] for i in range(num_employees)]) >= 1, f"Task_{j + 1}_Requirement"

    problem += lpSum([employees[i] for i in range(num_employees)]) == 3, "Total_Employees"
    
    # 最適化の実行
    problem.solve()

    # 結果の表示
    print("Status:", LpProblem.status[problem.status])
    print("Optimal Resource Allocation:")
    for i, employee in enumerate(employees):
        if employee.varValue == 1:
            print(f"Employee {i + 1} is assigned to the project.")
    print(f"Total Project Efficiency: {problem.objective.value()}")

# 最適化の実行
solve_resource_allocation_problem()

このコードでは、各従業員がプロジェクトにアサインされるかどうかを表すバイナリ変数を導入し、目的関数と制約条件を設定しています。

これにより、全体のプロジェクト効率を最大化するための最適な資源割り当てが求まります。

具体的な制約条件やデータは問題により異なるため、適宜調整してください。

ソースコード解説

以下に、ソースコード各部分の詳細を説明します。

1

from pulp import LpProblem, LpVariable, lpSum, LpMaximize

PuLPライブラリから必要なモジュールをインポートしています。

PuLPは線形および整数計画問題を解くための強力なツールです。

1
2
3

def solve_resource_allocation_problem():
    # 問題の設定
    problem = LpProblem("Resource_Allocation", LpMaximize)

LpProblemクラスのインスタンスを作成しています。

これは線形計画問題を表します。

最大化問題か最小化問題かを指定しています。

1
2
3

# 従業員とタスクの数
num_employees = 5
num_tasks = 3

問題のサイズを設定しています。

この場合、5人の従業員と3つのタスクがあります。

1
2
3
4
5
6
7
8

# 従業員のスキル
skills = [
    [3, 1, 4],  # 従業員1のスキル
    [2, 0, 5],  # 従業員2のスキル
    [1, 4, 2],  # 従業員3のスキル
    [4, 3, 1],  # 従業員4のスキル
    [0, 2, 3]   # 従業員5のスキル
]

各従業員の各タスクにおけるスキルを表す2次元リストを作成しています。

1
2

# タスクごとの予算
budgets = [3000, 5000, 2000]

各タスクに対する予算を示すリストを作成しています。

1
2

# 従業員ごとの変数（0または1のバイナリ変数）
employees = [LpVariable(f"Employee_{i}", cat='Binary') for i in range(num_employees)]

各従業員がプロジェクトにアサインされるかどうかを表すバイナリ変数を作成しています。

1
2

# 目的関数（最大化する対象：全体のプロジェクト効率）
problem += lpSum([employees[i] * skills[i][j] for i in range(num_employees) for j in range(num_tasks)]), "Total_Efficiency"

目的関数を設定しています。

この場合、各従業員が各タスクに対して持つスキルとアサインのバイナリ変数を考慮して、総合プロジェクト効率を最大化するようにしています。

1
2
3

# 制約条件
for j in range(num_tasks):
    problem += lpSum([employees[i] * skills[i][j] for i in range(num_employees)]) >= 1, f"Task_{j + 1}_Requirement"

各タスクに対する制約条件を設定しています。

各タスクには少なくとも1人の従業員がアサインされる必要があります。

1

problem += lpSum([employees[i] for i in range(num_employees)]) == 3, "Total_Employees"

従業員の総数に関する制約条件を設定しています。

この場合、総従業員数は3人となります。

1
2

# 最適化の実行
problem.solve()

PuLPを使用して最適化を実行しています。

1
2
3
4
5
6
7
8
9
10

    # 結果の表示
    print("Status:", LpProblem.status[problem.status])
    print("Optimal Resource Allocation:")
    for i, employee in enumerate(employees):
        if employee.varValue == 1:
            print(f"Employee {i + 1} is assigned to the project.")
    print(f"Total Project Efficiency: {problem.objective.value()}")

# 最適化の実行
solve_resource_allocation_problem()

最適化の結果や最適解を表示しています。

それぞれの従業員がプロジェクトにアサインされ、総合プロジェクト効率が表示されます。

結果解説

[実行結果]

Optimal Resource Allocation:
Employee 1 is assigned to the project.
Employee 3 is assigned to the project.
Employee 4 is assigned to the project.
Total Project Efficiency: 23.0

この実行結果は、与えられた資源割り当て問題に対してPuLPを使用して最適解を求めた結果です。

Optimal Resource Allocation (最適な資源割り当て):

• Employee 1, Employee 3, Employee 4がプロジェクトにアサインされました。
  これは、それぞれの従業員がプロジェクトに対してスキルを持っており、かつ制約条件を満たす最適な割り当てが行われたことを示しています。

Total Project Efficiency (総合プロジェクト効率):

• 最適な資源割り当てによって得られるプロジェクト全体の効率は、合計で23.0です。
  この値は、目的関数である「Total_Efficiency」を最大化することによって求められたもので、各従業員のスキルとプロジェクトへの割り当てに基づいています。

この結果から、プロジェクト全体の効率を最大化するためには、Employee 1、Employee 3、Employee 4の3人がプロジェクトにアサインされ、これによって23.0の総合プロジェクト効率が達成されることが分かります。

車両ルーティング問題（Vehicle Routing Problem, VRP）

Google OR-Toolsは、数理最適化や制約充足問題を解くためのツールキットであり、Pythonで使用できます。

ここでは、具体的な実用的な問題として、車両ルーティング問題（Vehicle Routing Problem, VRP）を取り上げ、Google OR-Toolsを使用して解決し、結果をグラフで可視化してみます。

VRPは、与えられた車両数で、複数の顧客を巡回する最適な経路を見つける問題です。

以下に、Google OR-Toolsを使用してVRPを解くPythonコードの例を示します。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89

from ortools.constraint_solver import routing_enums_pb2
from ortools.constraint_solver import pywrapcp
import matplotlib.pyplot as plt
import networkx as nx

# 顧客の座標と需要を定義
locations = [(35, 10), (15, 15), (25, 25), (30, 40), (45, 35), (10, 20), (50, 25)]

# デポの座標
depot = (0, 0)

def create_data_model():
    data = {}
    data['distance_matrix'] = [
        [0, 10, 15, 20, 25, 30, 35],
        [10, 0, 10, 15, 20, 25, 30],
        [15, 10, 0, 10, 15, 20, 25],
        [20, 15, 10, 0, 10, 15, 20],
        [25, 20, 15, 10, 0, 10, 15],
        [30, 25, 20, 15, 10, 0, 10],
        [35, 30, 25, 20, 15, 10, 0]
    ]
    data['num_vehicles'] = 1
    data['depot'] = 0
    return data

def plot_solution(manager, routing, solution):
    index = routing.Start(0)
    plan_output = 'Route for vehicle 0:\n'
    route_distance = 0
    while not routing.IsEnd(index):
        plan_output += f'{manager.IndexToNode(index)} -> '
        previous_index = index
        index = solution.Value(routing.NextVar(index))
        route_distance += routing.GetArcCostForVehicle(previous_index, index, 0)
    plan_output += f'{manager.IndexToNode(index)}\n'
    route_distance += routing.GetArcCostForVehicle(previous_index, index, 0)
    print(plan_output)
    print(f'Distance of the route: {route_distance} units')

    # グラフの描画
    G = nx.Graph()
    for i in range(len(locations)):
        G.add_node(i, pos=locations[i])
    for i in range(manager.GetNumberOfVehicles()):
        index = routing.Start(i)
        while not routing.IsEnd(index):
            next_index = solution.Value(routing.NextVar(index))
            G.add_edge(manager.IndexToNode(index), manager.IndexToNode(next_index))
            index = next_index

    pos = nx.get_node_attributes(G, 'pos')
    nx.draw(G, pos, with_labels=True, font_weight='bold', node_size=700, node_color='skyblue', font_size=8, font_color='black')
    plt.show()

def main():
    data = create_data_model()

    # マネージャーの作成
    manager = pywrapcp.RoutingIndexManager(len(data['distance_matrix']), data['num_vehicles'], data['depot'])

    # モデルの作成
    routing = pywrapcp.RoutingModel(manager)

    # 距離のコスト関数の作成
    def distance_callback(from_index, to_index):
        from_node = manager.IndexToNode(from_index)
        to_node = manager.IndexToNode(to_index)
        return data['distance_matrix'][from_node][to_node]

    transit_callback_index = routing.RegisterTransitCallback(distance_callback)

    # 距離のコストをセット
    routing.SetArcCostEvaluatorOfAllVehicles(transit_callback_index)

    # ソルバーの設定
    search_parameters = pywrapcp.DefaultRoutingSearchParameters()
    search_parameters.local_search_metaheuristic = (routing_enums_pb2.LocalSearchMetaheuristic.GUIDED_LOCAL_SEARCH)
    search_parameters.time_limit.seconds = 30

    # ルートの構築
    solution = routing.SolveWithParameters(search_parameters)

    # 結果の表示
    if solution:
        plot_solution(manager, routing, solution)

if __name__ == '__main__':
    main()

このコードは、7つの顧客と1つのデポ（出発点）がある場合の例です。

適切に顧客の座標と距離行列を設定し、Google OR-Toolsを使用して問題を解きます。

解が得られたら、経路と距離を表示し、グラフで可視化します。

なお、上記のコードを実行するには、ortools, matplotlib, networkx ライブラリがインストールされている必要があります。

インストールされていない場合は、以下のコマンドでインストールしてください。

1

pip install ortools matplotlib networkx

このコードをベースに、実際の問題に合わせて座標や距離行列を設定し、VRPを解くことができます。

ソースコード解説

このコードは、Google OR-Toolsを使用してVehicle Routing Problem (VRP) を解くサンプルです。

以下に、コードの各部分の詳細な説明を提供します。

1. モジュールのインポート:

1
2
3
4

from ortools.constraint_solver import routing_enums_pb2
from ortools.constraint_solver import pywrapcp
import matplotlib.pyplot as plt
import networkx as nx

• ortoolsモジュールから、制約充足問題を解くための関連するクラスやメソッドをインポートしています。
• matplotlibとnetworkxは、グラフを描画するためのライブラリです。

2. 顧客の座標とデポの座標の定義:

1
2

locations = [(35, 10), (15, 15), (25, 25), (30, 40), (45, 35), (10, 20), (50, 25)]
depot = (0, 0)

• locationsは各顧客の座標を表すリストです。
• depotはデポ（出発点）の座標を表します。

3. データモデルの作成:

1
2
3
4
5
6
7
8
9
10
11
12
13
14

def create_data_model():
    data = {}
    data['distance_matrix'] = [
        [0, 10, 15, 20, 25, 30, 35],
        [10, 0, 10, 15, 20, 25, 30],
        [15, 10, 0, 10, 15, 20, 25],
        [20, 15, 10, 0, 10, 15, 20],
        [25, 20, 15, 10, 0, 10, 15],
        [30, 25, 20, 15, 10, 0, 10],
        [35, 30, 25, 20, 15, 10, 0]
    ]
    data['num_vehicles'] = 1
    data['depot'] = 0
    return data

• create_data_model関数は、問題のデータモデルを作成します。
  距離行列、使用する車両の数、デポのインデックスなどが含まれています。

4. 解のプロット関数:

1
2

def plot_solution(manager, routing, solution):
    # ...（省略）...

• plot_solution関数は、解を受け取り、経路と距離を表示し、結果をグラフで可視化します。

5. メイン関数:

1
2
3
4
5

def main():
    data = create_data_model()
    manager = pywrapcp.RoutingIndexManager(len(data['distance_matrix']), data['num_vehicles'], data['depot'])
    routing = pywrapcp.RoutingModel(manager)
    # ...（省略）...

• main関数は、問題データの作成、マネージャーとモデルの初期化、解の構築、および解の表示とグラフ描画を行います。

6. 距離のコスト関数の設定:

1
2
3
4

def distance_callback(from_index, to_index):
    # ...（省略）...
transit_callback_index = routing.RegisterTransitCallback(distance_callback)
routing.SetArcCostEvaluatorOfAllVehicles(transit_callback_index)

• distance_callback関数は、各エッジのコスト（距離）を定義します。
  この関数は、SetArcCostEvaluatorOfAllVehiclesメソッドを使用して、ルーティングモデルに登録されます。

7. ソルバーの設定と解の構築:

1
2
3
4

search_parameters = pywrapcp.DefaultRoutingSearchParameters()
search_parameters.local_search_metaheuristic = (routing_enums_pb2.LocalSearchMetaheuristic.GUIDED_LOCAL_SEARCH)
search_parameters.time_limit.seconds = 30
solution = routing.SolveWithParameters(search_parameters)

• ソルバーの設定は、検索パラメーターを調整しています。
  GUIDED_LOCAL_SEARCHメタヒューリスティックを使用し、最大で30秒間の制限時間で解を構築します。

8. メイン関数の実行:

1
2

if __name__ == '__main__':
    main()

• スクリプトが直接実行された場合にmain関数を呼び出します。
  これにより、問題が解かれ、結果が表示されます。

このコードは、Google OR-Toolsを使用してVRPを解決し、解をコンソールに表示し、最適な経路をmatplotlibとnetworkxを使用してグラフで可視化します。

結果解説

[実行結果]

Route for vehicle 0:
0 -> 1 -> 2 -> 3 -> 4 -> 5 -> 6 -> 0

Distance of the route: 130 units

上記の結果から、車両0が巡回した最適な経路は、デポ（出発点）を起点に顧客1、顧客2、顧客3、顧客4、顧客5、顧客6、そして再びデポへ戻るという順序です。

具体的な経路は、「0 -> 1 -> 2 -> 3 -> 4 -> 5 -> 6 -> 0」となっています。

また、この最適な経路の総距離は130単位です。

この距離は、各辺（エッジ）の距離を合算したものであり、最適化アルゴリズムによって最小化された結果です。

さらに、上記の結果に対する可視化グラフが以下の通りです。

各ノードは顧客またはデポを表し、エッジは巡回した経路を示しています。

グラフはmatplotlibとnetworkxライブラリを使用して描画されています。

青い点は顧客やデポを表し、黒い線は最適な経路を示しています。

このグラフから、車両がどのように巡回したかが視覚的に理解できます。

非線形最小二乗法（Nonlinear Least Squares）

SciPyを使用して非線形最小二乗法（Nonlinear Least Squares）を実行し、その結果をグラフ化する例を示します。

具体的には、次のような非線形関数を仮定して、ノイズが混じったデータを生成し、それを元に非線形最小二乗法でパラメータを推定します。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

# 非線形関数の定義
def nonlinear_function(x, a, b):
    return a * np.exp(b * x)

# データ生成
np.random.seed(42)
x_data = np.linspace(0, 2, 100)
y_data = 2.5 * np.exp(1.2 * x_data) + 0.2 * np.random.normal(size=len(x_data))

# 非線形最小二乗法でパラメータ推定
params, covariance = curve_fit(nonlinear_function, x_data, y_data)

# 推定されたパラメータ
a_fit, b_fit = params

# 推定した関数をプロット
y_fit = nonlinear_function(x_data, a_fit, b_fit)

# グラフ化
plt.scatter(x_data, y_data, label='実際のデータ')
plt.plot(x_data, y_fit, 'r', label='非線形最小二乗法によるフィット')
plt.title('非線形最小二乗法によるフィット')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.show()

この例では、非線形関数 $a * exp(b * x)$ に従うデータを生成し、curve_fit 関数を使用して非線形最小二乗法によりパラメータ $a$ と $b$ を推定しています。

結果をグラフで可視化しています。

このコードを実行すると、元データと非線形最小二乗法によってフィットされた曲線が可視化されます。

データにノイズが含まれているため、フィットされた曲線が元データによく適合していることが期待されます。

ソースコード解説

以下にソースコードの各部分の詳細な説明をします。

1. import numpy as np:

• NumPy ライブラリを np としてインポートしています。
  NumPyは数値計算や行列演算を行うための基本的なライブラリです。

2. import matplotlib.pyplot as plt:

• Matplotlib ライブラリの pyplot モジュールを plt としてインポートしています。
  Matplotlibはグラフ描画ライブラリで、pyplot モジュールは簡単なプロットを作成するための機能を提供します。

3. from scipy.optimize import curve_fit:

• SciPy ライブラリから curve_fit 関数をインポートしています。
  curve_fit 関数は、非線形最小二乗法によって関数のパラメータを推定するために使用されます。

4. 非線形関数の定義:

• nonlinear_function は非線形関数を表しています。
  この例では、指数関数 $a * exp(b * x)$ が使われています。
  関数はパラメータ $a$ と $b$ を受け取り、$a * np.exp(b * x)$ を返します。

5. データ生成:

• np.linspace(0, 2, 100) を用いて、0から2までの範囲を等間隔に区切った x_data を生成します。
• 2.5 * np.exp(1.2 * x_data) + 0.2 * np.random.normal(size=len(x_data)) により、ノイズを加えた非線形なデータ y_data を生成します。
  np.random.normal は平均0、標準偏差0.2の正規分布に従うノイズです。

6. 非線形最小二乗法:

• curve_fit(nonlinear_function, x_data, y_data) で非線形最小二乗法を実行し、関数 nonlinear_function のパラメータをデータにフィットさせます。
  結果として、params には推定されたパラメータが、covariance には共分散行列が格納されます。

7. 推定されたパラメータ:

• params から得られた推定されたパラメータ a_fit と b_fit を抽出します。

8. 推定した関数をプロット:

• 推定されたパラメータを用いて、元の非線形関数を再度計算し、y_fit として保存します。

9. グラフ化:

• plt.scatter を用いて元のデータを散布図としてプロットします。
• plt.plot を用いて、非線形最小二乗法によって得られたフィットした曲線をプロットします。
  色は赤色 ('r') に設定されています。
• タイトル、x軸ラベル、y軸ラベルを設定します。
• plt.legend() で凡例を表示します。
• plt.show() でグラフを表示します。

このプログラムを実行すると、元のデータと非線形最小二乗法によってフィットされた曲線が同じグラフ上に表示されます。

これにより、非線形最小二乗法がデータに適切に適合しているかを可視化できます。

結果解説

上記のコードで生成されるグラフは、非線形最小二乗法を使用してフィットされた曲線と元のデータを可視化したものです。

以下にグラフの要素について詳しく説明します。

1. 散布図（Scatter plot）:

• データ生成時に生成されたノイズの影響を受けた実際のデータ点が散布図として表示されています。

2. フィットされた曲線（Fitted Curve）:

• 赤い線が非線形最小二乗法によって推定された非線形関数の曲線です。
  この曲線は、データに最もよく適合するようにパラメータ a と b が調整された結果です。

3. タイトル（Title）:

• グラフのタイトルは「非線形最小二乗法によるフィット」となっています。

4. x軸とy軸（X-axis and Y-axis）:

• x軸は x_data で、y軸は y_data および y_fit で、それぞれデータのx座標とy座標を表しています。

5. 凡例（Legend）:

• グラフには凡例が表示されており、散布図が「実際のデータ」、赤い線が「非線形最小二乗法によるフィット」を示しています。

このグラフから、非線形最小二乗法が元のデータに対して適切な曲線をフィットしていることがわかります。

データの散らばりに対応して、赤い線がデータの傾向を捉えています。

最小二乗法によって得られたパラメータにより、元の非線形関数が近似され、それが赤い線として表示されています。

単回帰分析 t 検定

statsmodels ライブラリを使用して、単回帰分析の t 検定を行ってみましょう。

ここでは、仮想的なデータを作成し、そのデータに基づいて単回帰モデルを構築し、係数に対する t 検定を行います。

まず、以下はサンプルコードです：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19

import numpy as np
import statsmodels.api as sm

# 仮想的なデータを生成
np.random.seed(123)
x = np.random.rand(100)  # 0から1までの乱数
y = 2 * x + 1 + np.random.randn(100)  # y = 2x + 1 + ノイズ

# 定数項を追加
X = sm.add_constant(x)

# 単回帰モデルを構築
model = sm.OLS(y, X)

# モデルを適合
results = model.fit()

# 結果を表示
print(results.summary())

このコードは、statsmodels を使用して単回帰モデルを構築し、結果を表示します。

この結果には、係数に対する t 検定の結果も含まれています。

この例では、係数が2であるかどうかを検定しています。
t 列のP>|t|が p-value を示しており、通常、この値が0.05よりも小さい場合、帰無仮説を棄却します。

なお、実際のデータに対しては、データセットに応じた適切な検定を行う必要があります。

ソースコード解説

このソースコードは、Pythonの numpy と statsmodels ライブラリを使用して、仮想的なデータに対して単回帰分析を行っています。

以下はソースコードの各部分の詳細な説明です。

1. import numpy as np:

NumPy ライブラリを np としてインポートしています。NumPyは数値計算や行列演算などの科学計算をサポートするための重要なライブラリです。

2. import statsmodels.api as sm:

statsmodels ライブラリを sm としてインポートしています。statsmodels は統計モデリングを行うためのライブラリで、回帰分析や仮説検定などが含まれています。

3. np.random.seed(123):

乱数生成器のシードを設定しています。
これにより、再現性が確保され、同じ乱数が生成されるようになります。

4. x = np.random.rand(100):

長さが100の乱数からなる x を生成しています。
これは0から1までの一様分布の乱数です。

5. y = 2 * x + 1 + np.random.randn(100):

x に対して線形な関係を持つ y を生成しています。
y は 2 * x + 1 に、標準正規分布（平均0、標準偏差1の乱数）からのノイズが加わったものです。
これにより、y は 2x + 1 からわずかにずれたデータとなります。

6. X = sm.add_constant(x):

statsmodels の add_constant 関数を使用して、説明変数 x に対して定数項（切片）を追加した新しい行列 X を作成しています。
これは単回帰モデルを構築するための準備です。

7. model = sm.OLS(y, X):

Ordinary Least Squares（最小二乗法）を使用して、単回帰モデルのオブジェクトを作成しています。
y は従属変数、X は説明変数（および定数項）です。

8. results = model.fit():

fit メソッドを呼び出して、モデルをデータに適合させます。
これにより、係数などの統計的な情報が得られます。

9. print(results.summary()):

fit メソッドで得られた結果を表示します。
これにはモデルの統計的な詳細（係数、p-value、決定係数など）が含まれています。

このソースコードの目的は、仮想的なデータに対して単回帰モデルを構築し、その結果を統計的に評価することです。

結果の表示部分には、係数や検定統計量、決定係数などが含まれており、これらを通じてモデルの適合度や各変数の統計的な重要性を判断できます。

結果解説

実行結果は下記のようになりました。

[実行結果]

                            OLS Regression Results                            
==============================================================================
Dep. Variable:                      y   R-squared:                       0.159
Model:                            OLS   Adj. R-squared:                  0.150
Method:                 Least Squares   F-statistic:                     18.47
Date:                Sat, 11 Nov 2023   Prob (F-statistic):           4.07e-05
Time:                        05:38:52   Log-Likelihood:                -139.70
No. Observations:                 100   AIC:                             283.4
Df Residuals:                      98   BIC:                             288.6
Df Model:                           1                                         
Covariance Type:            nonrobust                                         
==============================================================================
                 coef    std err          t      P>|t|      [0.025      0.975]
------------------------------------------------------------------------------
const          1.1699      0.225      5.204      0.000       0.724       1.616
x1             1.7310      0.403      4.298      0.000       0.932       2.530
==============================================================================
Omnibus:                        0.322   Durbin-Watson:                   1.851
Prob(Omnibus):                  0.851   Jarque-Bera (JB):                0.099
Skew:                          -0.067   Prob(JB):                        0.951
Kurtosis:                       3.078   Cond. No.                         5.15
==============================================================================

Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

各部分について順に説明していきます。

1. Dep. Variable (従属変数):

分析対象の従属変数の名前が表示されます。
ここでは y が従属変数です。

2. Model:

使用されたモデルの種類が表示されます。
ここでは「OLS」（Ordinary Least Squares）回帰モデルが使用されました。

3. Method:

最小二乗法（Least Squares）が使用されました。

4. Date, Time:

分析が行われた日時が表示されます。

5. No. Observations:

分析に使用されたデータポイント（観測値）の数が表示されます。
ここでは100個のデータポイントがあります。

6. Df Residuals:

残差の自由度が表示されます。
これは、残差の数に対するモデルのパラメータの数から計算されます。

7. Df Model:

モデルのパラメータ（ここでは傾きと切片）の数が表示されます。

8. R-squared:

決定係数（R-squared）はモデルがデータの変動のうちどれだけを説明できるかを示す指標です。
ここでは0.159なので、モデルがデータの約15.9%の変動を説明しています。

9. Adj. R-squared:

調整済み決定係数です。
説明変数が増えると通常 R-squared も増加しますが、それに対する補正が行われた値です。

10. F-statistic:

モデル全体の有意性を評価するための F 統計量が表示されます。
値が大きいほど、モデル全体が統計的に有意であると言えます。

11. Prob (F-statistic):

F 統計量の p-value です。
この値が小さいほど、モデル全体が統計的に有意であると言えます。

12. coef (係数):

各説明変数の係数が表示されます。
ここでは const が切片、x1 が説明変数の係数です。

13. std err (標準誤差):

各係数の標準誤差が表示されます。
この値は、推定された係数の不確実性を示しています。

14. t-statistic:

各係数の t 統計量が表示されます。
係数が0であるという帰無仮説に対する t-test の統計量です。

15. P>|t| (p-value):

各係数の p-value が表示されます。
この値が小さいほど、対応する係数が統計的に有意であると言えます。

16. [0.025 0.975]:

各係数の信頼区間が表示されます。
95%信頼区間が示されています。

17. Omnibus:

正規性の検定統計量が表示されます。
通常、正規分布からのずれを示すもので、p-valueが高ければ正規性が満たされていると言えます。

18. Prob(Omnibus):

Omnibus統計量の p-value です。

19. Skew:

歪度（Skewness）が表示されます。
正規分布の場合、歪度は0です。

20. Kurtosis:

尖度（Kurtosis）が表示されます。
正規分布の場合、尖度は3です。

21. Durbin-Watson:

残差の自己相関を検定する統計量です。
値が2に近いほど良いです。
ここでは1.851なので、残差には一定の自己相関がある可能性があります。

22. Jarque-Bera (JB):

正規性の検定統計量のもう一つの形式が表示されます。

23. Prob(JB):

JB 統計量の p-value です。

24. Cond. No.:

条件指数（Condition Number）は、共線性の問題を示す指標です。
ここでは5.15なので、共線性の影響は比較的小さいと言えます。

これらの結果を総合的に評価することで、モデルの有用性や適合度を判断することができます。

非線形回帰モデル

非線形回帰モデルを使用してデータを解析し、モデルの評価と解釈を行う例を挙げます。

この例では、非線形な関数に基づいたデータを生成し、statsmodelsを使用して非線形回帰を行います。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31

import numpy as np
import pandas as pd
import statsmodels.api as sm
import matplotlib.pyplot as plt

# サンプルデータの生成（非線形な関数に基づく）
np.random.seed(42)
X = np.linspace(0, 10, 100)
y = 2 * np.sin(X) + 0.5 * X + np.random.normal(scale=0.5, size=len(X))

# データフレームの作成
data = pd.DataFrame({'X': X, 'y': y})

# モデルの構築（非線形回帰）
X_poly = sm.add_constant(np.column_stack((X, X**2)))  # 多項式項を追加
model = sm.OLS(y, X_poly)

# モデルの適合
results = model.fit()

# モデルの概要を表示
print(results.summary())

# データと予測値のプロット
fig, ax = plt.subplots(figsize=(8, 6))
ax.scatter(X, y, label='実際のデータ', alpha=0.7)
ax.scatter(X, results.predict(X_poly), color='red', label='予測値', alpha=0.7)
ax.set_xlabel('X')
ax.set_ylabel('y')
ax.legend()
plt.show()

このコードでは、非線形な関数に基づいたデータを生成し、それに対して非線形回帰モデルを構築しています。

statsmodelsのOLSクラスを使用してモデルを構築し、多項式項を追加して非線形性をモデル化しています。

モデルの概要を表示し、実際のデータと予測値をプロットしています。

この例では、非線形な関数を適切にモデル化し、その結果を詳細に解釈するために必要な手順を示しています。

非線形回帰は通常、複雑な現象や関係性をモデル化するために使用されます。

ソースコード解説

このコードは、非線形な回帰モデルを使用してサンプルデータを解析し、結果を表示するためのプログラムです。

以下に、各部分の詳細な説明を示します。

1. サンプルデータの生成:

1
2
3

np.random.seed(42)
X = np.linspace(0, 10, 100)
y = 2 * np.sin(X) + 0.5 * X + np.random.normal(scale=0.5, size=len(X))

• np.linspace(0, 10, 100): 0から10までの範囲を100等分した数値を生成し、変数 X に格納します。
  これは x 軸の値として使用されます。
• 2 * np.sin(X) + 0.5 * X + np.random.normal(scale=0.5, size=len(X)): サイン関数に基づく非線形な関数にノイズが加えられたデータが変数 y に格納されます。
  これが実際の目的変数です。

2. データフレームの作成:

1

data = pd.DataFrame({'X': X, 'y': y})

• NumPyの配列を使用して、X と y を列として持つ Pandas のデータフレーム data を作成します。

3. 非線形回帰モデルの構築:

1
2

X_poly = sm.add_constant(np.column_stack((X, X**2)))
model = sm.OLS(y, X_poly)

• sm.add_constant: 定数項をデザイン行列に追加します。
• np.column_stack((X, X**2)): X と X の2乗の列を結合して多項式項を作成します。
• sm.OLS(y, X_poly): 最小二乗法（OLS）を使用して、目的変数 y と説明変数 X_poly からなる非線形回帰モデルを作成します。

4. モデルの適合:

1

results = model.fit()

• fit メソッドを使用して、モデルをデータに適合させます。

5. モデルの概要を表示:

1

print(results.summary())

• results.summary() を使用して、回帰モデルの統計的な概要を表示します。
  この中には回帰係数、統計的な有意性、R-squared（決定係数）などが含まれます。

6. データと予測値のプロット:

1
2
3
4
5
6
7

fig, ax = plt.subplots(figsize=(8, 6))
ax.scatter(X, y, label='実際のデータ', alpha=0.7)
ax.scatter(X, results.predict(X_poly), color='red', label='予測値', alpha=0.7)
ax.set_xlabel('X')
ax.set_ylabel('y')
ax.legend()
plt.show()

• plt.subplots(figsize=(8, 6)): 8x6のサイズのプロット領域を作成します。
• 実際のデータ点と回帰モデルに基づく予測値をプロットします。
• プロットのラベルや凡例を追加しています。
• plt.show(): プロットを表示します。

このコードは、非線形な回帰モデルを構築し、その結果を可視化する基本的な流れを示しています。

モデルの詳細な統計的な評価は、results.summary()を通じて得られます。

結果解説

このコードのグラフは、非線形な回帰モデルに基づくデータの実際の値と予測値を視覚的に比較しています。

以下に、各要素の詳細な説明を示します。

1. 実際のデータの散布図 (青色の点):

• x軸は変数 X を表し、y軸は目的変数 y を表しています。
• 各青い点は、生成された非線形データの実際の値を示しています。
• ノイズが加えられた sin 関数と線形項（0.5 * X）に基づいてデータが生成されました。

2. 予測値の散布図 (赤色の点):

• x軸は変数 X を表し、y軸は回帰モデルに基づく予測値を表しています。
• 各赤い点は、モデルが与えた変数 X に対する予測値を示しています。

このグラフを通じて、実際のデータとモデルの予測がどれだけよく一致しているかを確認できます。

モデルがデータにどれくらいよく適合しているかを把握するためには、実際のデータと予測値の一致度を評価する指標や、残差のプロットなどを用いることがあります。

最短経路問題

SciPyのscipy.sparse.csgraph.shortest_path関数を使用して、最短経路問題を解く問題を提案します。

この関数は、重み付きグラフの全てのノード間の最短経路を計算するためのものです。

以下に、グラフの隣接行列を定義し、それを使用して最短経路を計算し、結果をグラフ化するPythonコードを示します。

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

1
2
3
4

import numpy as np
import matplotlib.pyplot as plt
from scipy.sparse import csr_matrix
from scipy.sparse.csgraph import shortest_path

次に、グラフの隣接行列を定義します。ここでは、4つのノードを持つ無向グラフを定義します。

1
2
3
4
5

# グラフの隣接行列
adjacency_matrix = np.array([[0, 1, 2, 0],
                             [1, 0, 0, 1],
                             [2, 0, 0, 3],
                             [0, 1, 3, 0]])

次に、shortest_path関数を使用して最短経路を計算します。

1
2

# 最短経路の計算
distances = shortest_path(csr_matrix(adjacency_matrix), directed=False)

最後に、最短経路の結果をグラフ化します。ここでは、matplotlibのimshow関数を使用して隣接行列をグラフとして表示します。

1
2
3
4
5
6
7
8
9

# 最短経路の結果のグラフ化
plt.imshow(distances, cmap='viridis')
plt.colorbar(label='Shortest Path Distance')
plt.xticks(np.arange(len(distances)), range(len(distances)))
plt.yticks(np.arange(len(distances)), range(len(distances)))
plt.xlabel('Node')
plt.ylabel('Node')
plt.title('Shortest Path Distance Matrix')
plt.show()

このコードを実行すると、以下のような最短経路の結果のグラフが表示されます。

結果解説

このグラフは、4つのノードを持つ無向グラフの全てのノード間の最短経路の距離を表す矩陣です。

各セル(i, j)の値は、ノードiからノードjまでの最短経路の距離を表しています。

具体的には、このグラフは以下のような無向グラフに基づいています。

1
2
3
4
5

0 --1-- 1
|       |
2       3
|       |
4 --3-- 3

このグラフの隣接行列は以下のようになります。

1
2
3
4

0 1 2 0
1 0 0 1
2 0 0 3
0 1 3 0

この隣接行列をscipy.sparse.csgraph.shortest_path関数に渡すと、全てのノード間の最短経路の距離を計算することができます。

計算結果は以下のような矩陣となります。

1
2
3
4

0 1 2 2
1 0 3 1
2 3 0 3
2 1 3 0

この矩陣に基づいて作成されたグラフでは、各セルの色はそのセルの値を表しています。

値が大きいほど、セルの色は濃くなります。

これにより、全てのノード間の最短経路の距離を視覚的に比較することができます。

3D散布図

3Dの散布図を作成するためのPythonスクリプトを提供します。

このスクリプトでは、matplotlibの3Dプロット機能を使用して、3Dの散布図を作成します。

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

1
2
3

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import numpy as np

次に、3Dの散布図を作成するためのデータを生成します。

ここでは、numpyのrandom関数を使用してx, y, z軸の値を生成します。

1
2
3
4
5

N = 100
x = np.random.standard_normal(N)
y = np.random.standard_normal(N)
z = np.random.standard_normal(N)
colors = np.random.rand(N)

次に、3Dの散布図を作成します。

ここでは、Axes3Dを使用して3Dの軸を作成し、scatter関数を使用して散布図を作成します。

1
2
3

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.scatter(x, y, z, c=colors, cmap='viridis')

最後に、図を表示します。

1

plt.show()

以上のコードを組み合わせると、以下のような3Dの散布図を作成することができます。