ChatGPTとPythonで学ぶ 非尺度型多次元尺度構成法 (NMDS)



皆さんこんにちは。今日も引き続きChatGPT先生をお迎えして、「ChatGPTとPythonで学ぶ 非尺度型多次元尺度構成法」というテーマで雑談したいと思います。それではChatGPT先生、よろしくお願いします。

---

### **assistant**:

こんにちは！引き続きよろしくお願いします。非尺度型多次元尺度構成法について学びましょう。

まず、非尺度型多次元尺度構成法（Non-metric Multidimensional Scaling、NMDSと略すこともあります）を簡単に説明します。これは、個体間の相対的な距離を保存しながら高次元データを低次元に圧縮する手法です。つまり、距離の大小関係だけを保持して、個体間の絶対的な距離は無視します。これにより、非メトリック空間でのデータのパターンを調査することが可能となります。

それでは、Pythonを使って非尺度型多次元尺度構成法を実行してみましょう。Pythonでは、主に `sklearn` パッケージの `manifold` モジュールにある `MDS` クラスを使用します。

※注意点としては、NMDSは計算量が大きいため、データ数が大きい場合は計算に時間がかかる可能性があります。また、初期値による影響が大きいため、結果は繰り返し実行するたびに変わることがあります。

---

### **user**:

MDSとNMDSの違いは何でしょうか？詳しく教えてください。

---

### **assistant**:

よい質問です！では、MDS (Metric Multidimensional Scaling)とNMDS (Non-metric Multidimensional Scaling) の違いについて説明します。

MDSとNMDSはどちらも多次元スケーリングの手法ですが、対象となる距離データの種類と解釈の方法に違いがあります。

1. **MDS** (Metric Multidimensional Scaling):
MDSはユークリッド距離やマンハッタン距離などの「メトリック距離」に基づいて次元削減を行います。MDSでは、元の高次元空間での距離と、削減後の低次元空間での距離の差（ストレス）を最小化するように次元削減を行います。そのため、MDSによって得られるプロットでは、オブジェクト間の距離ができる限り忠実に反映されます。

2. **NMDS** (Non-metric Multidimensional Scaling)
一方、NMDSは任意の相互関係（メトリックではない距離）を扱うことができます。NMDSでは、高次元空間と低次元空間でオブジェクト間の順序尺度レベルの関係性（つまり距離のランキング）が保存されるように次元削減を行います。つまり、個体間の絶対的な距離よりも、個体間の距離の順番（ランキング）を保持することに重点を置きます。

これらの違いから、MDSとNMDSはそれぞれ異なる種類のデータと問いに対して適用されます。たとえばMDSは主に量的な相互作用を観察するために使用され、個体間の絶対的な距離を解釈するのに有用です。一方、NMDSは質的または半量的な相互作用を観察するために使用され、個体間の距離関係の順位（ランキング）を解釈するのに有用です。

---

### **user**:

それでは、NMDSについて学ぶためのPythonコードを書いてください。何かのデータに対してNMDSを実行して結果を可視化してください。

---

### **assistant**:

はい、それではPythonで非尺度型多次元尺度構成法を実行するサンプルコードをご紹介します。以下に示すコードでは、アイリスのデータセットに対して非尺度型多次元尺度構成法を用いてデータを2次元に削減し、その結果をプロットします。

なお、PythonのsklearnライブラリではNMDSの実装が直接提供されていないため、`MDS`クラスに`dissimilarity='precomputed'`オプションを用いて非尺度型多次元尺度構成法をエミュレートします。これはオブジェクト間の距離行列を直接入力とするためです。

```python
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.metrics import pairwise_distances
from sklearn.manifold import MDS

# Load Iris dataset
iris = datasets.load_iris()
X = iris.data
y = iris.target

# Compute pairwise distance matrix
D = pairwise_distances(X)

# Apply non-metric MDS
model = MDS(n_components=2, dissimilarity='precomputed', random_state=1)
X_transformed = model.fit_transform(D)

# Plot results
plt.scatter(X_transformed[:, 0], X_transformed[:, 1], c=y)
plt.show()
```

![output_0_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/211162/946e75c3-2c2b-24cb-b929-ccfbee929d63.png)


このコードでは、まずアイリスのデータセットをロードします。そのデータセットからペアワイズの距離行列を計算し、非尺度型多次元尺度構成法を適用します。そして、変換された2次元データをプロットします。プロットの結果は、3種類のアイリスの花（yの値で色分け）が2次元空間にどのように配置されているかを示します。

なお、このコードは似たオブジェクトが高次元空間で近くなることを目指しています。しかし、尺度型のMDSと異なり、非尺度型のMDSは原空間の距離の順位（相対的な距離の大小）を維持することを目指します。これにより、相対的な距離が一定でも、原空間での絶対的な距離に基づくバイアスの影響を受けにくくなります。

---

### **user**:

pairwise_distances って何でしたっけ

---

### **assistant**:

`pairwise_distances`は、Scikit-learn（`sklearn`）ライブラリに含まれる関数です。この関数は、データセット中の全てのペアに対して指定した距離計算を行います。

具体的にの使い方は以下の通りです。

```python
from sklearn.metrics import pairwise_distances

# Example data
X = [[1, 2], [3, 4], [5, 6]]

# Compute pairwise distances
D = pairwise_distances(X)
```

上記の例では、`D`は`X`の全ての行間のペアワイズ距離を含む行列になります。つまり、`D[i, j]`は`X[i]`と`X[j]`の距離となります。

距離計算の方法（メトリクス）はデフォルトではユークリッド距離が用いられますが、他の距離尺度（例えばマンハッタン距離、コサイン距離など）を用いることもできます。これは、`pairwise_distances`関数の`metric`引数を用いて指定します。

非尺度型多次元尺度構成法（NMDS）を実装する際には、`pairwise_distances`関数で元のデータのペアワイズ距離を計算し、その距離行列を`MDS`クラスの`fit_transform`メソッドに渡します。

以上が`pairwise_distances`の基本的な説明となります。何か他にご質問がありましたらお知らせください。

---

### **user**:

MDSとNMDSの違いは、距離の定義だけですか？

---

### **assistant**:

私が以前言及したように、MDSとNMDSの主な違いは、それらが処理する距離データとその解釈の方向性にあります。しかし、その背後にある数学的な最適化問題（すなわち、損失関数や目的関数）にも違いがあります。

- **MDS（Metric Multidimensional Scaling）:** MDSは、多次元空間から低次元空間へのマッピングが原始的な距離をできるだけ維持するようにします。つまり、MDSは低次元空間での距離と元の距離との間の差（ストレス）を最小化しようとします。それは「尺度レベル」の情報を保持しようとするため、これを「メトリック」または「尺度型」の方法と呼びます。

- **NMDS（Non-metric Multidimensional Scaling）:** NMDSの目標は、距離の「順序」を維持することです。つまり、オブジェクト間の相対距離の順序（すなわち、「AはBよりもCに近い」など）が低次元空間でも維持されることを最適化します。そのため、NMDSは順序レベルの情報を保持しようとするため、これを「非メトリック」または「非尺度型」の方法と呼びます。

以上がMDSとNMDSの主な違いです。そのため、適用される問題や利用されるデータの性質によって、MDSまたはNMDSのどちらを選択するかが決まります。例えば、距離が順序情報を重視する非尺度データ（例: 初級、中級、上級などのカテゴリ）の場合、NMDSの方が適切です。

---

### **user**:

イリスのデータに対して、MDSとNMDSを実行し、結果を可視化して比較できるようなコードを書いてください。

---

### **assistant**:

イリスデータセットに対してまずMDS（Metric Multidimensional Scaling）を適用し、次にNMDS（Non-Metric Multidimensional Scaling）を適用して、その結果を比較するためのPythonコードを以下に示します。結果は2次元の散布図として表示します。

```python
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.manifold import MDS
from sklearn.metrics import pairwise_distances

# Load Iris dataset
iris = datasets.load_iris()
X = iris.data
y = iris.target

# Apply MDS
mds = MDS(n_components=2, random_state=1)
X_mds = mds.fit_transform(X)

# Compute pairwise distance matrix
D = pairwise_distances(X)

# Apply non-metric MDS
nmds = MDS(n_components=2, dissimilarity='precomputed', random_state=1)
X_nmds = nmds.fit_transform(D)

# Plot results
fig, ax = plt.subplots(1, 2, figsize=(12, 6))

# MDS plot
ax[0].scatter(X_mds[:, 0], X_mds[:, 1], c=y)
ax[0].set_title('MDS')

# NMDS plot
ax[1].scatter(X_nmds[:, 0], X_nmds[:, 1], c=y)
ax[1].set_title('NMDS')

plt.show()
```

![output_2_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/211162/d07172b0-346d-3ef5-c5cd-68fafdc3457b.png)


このコードでは、MDSとNMDSの結果を並べてプロットします。プロットの色は各アイリスの種類によって異なります。これにより、これら2つの方法が同じデータセットにどのように影響を与えるかを視覚的に比較することができます。

ただし、再度強調しますが、MDSとNMDSの結果は初期値に依存するため、計算ごとに結果が変わる可能性があります。

なお、`random_state`パラメータを設定することで、実行するたびに同じ結果が得られるようにすることができます。これは結果の再現性を確保するために重要なステップです。

---

### **user**:

うーん、MDSとNMDSで結果が違うようには見えないですね。

---

### **assistant**:

ご観察の通り、イリスのデータセットを用いたこの特定の例では、MDSとNMDSの結果は視覚的には類似しています。実際、あるデータセットに対しては、MDSとNMDSが非常に類似した結果をもたらすことがあります。

ただし、MDSとNMDSが同じように見える場合でも、それらが目指す目的は異なります。MDSはオブジェクト間の距離を可能な限り忠実に再現しようとします。一方、NMDSはオブジェクト間の距離のランキング（順序）を可能な限り保持しようとします。

したがって、この2つの手法は本質的には異なる問題を解決しようとしています。これらの違いが重要となるのは、距離尺度（MDS）と距離の順序（NMDS）の両方がデータの解釈に大きな影響を与えるようなケースです。

このような状況では、求められる情報の種類（すなわち、尺度情報または順序情報）に基づいて、MDSまたはNMDSを選択することになるでしょう。

---

### **user**:

ペンギンのデータセットに対して、MDSとNMDSを実行し、結果を可視化して比較できるようなコードを書いてください。

---

### **assistant**:

はい、そのコードを準備します。ただし、まず`seaborn`ライブラリを利用してペンギンのデータセットを読み込みます。`seaborn`ライブラリがインストールされていない場合は、`!pip install seaborn`を実行してインストールしてください。また、ペンギンのデータセットは数値データとカテゴリデータの両方を含むため、数値データのみを取り出します。

以下のコードでは、ペンギンデータセットの数値データに対してまずMDSを適用し、次にNMDSを適用し、それぞれの結果をプロットします。

```python
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.preprocessing import StandardScaler
from sklearn.manifold import MDS
from sklearn.metrics import pairwise_distances

# Load Penguin dataset
penguins = sns.load_dataset('penguins').dropna().reset_index(drop=True)
X = penguins.select_dtypes(include=['float64', 'int64'])
y = penguins['species']

# Standardize the feature vector
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)

# Apply MDS
mds = MDS(n_components=2, random_state=1)
X_mds = mds.fit_transform(X_scaled)

# Compute pairwise distance matrix
D = pairwise_distances(X_scaled)

# Apply non-metric MDS
nmds = MDS(n_components=2, dissimilarity='precomputed', random_state=1)
X_nmds = nmds.fit_transform(D)

# Plot results
fig, ax = plt.subplots(1, 2, figsize=(12, 6))

# MDS plot
scatter = ax[0].scatter(X_mds[:, 0], X_mds[:, 1], c=y.astype('category').cat.codes)
ax[0].set_title('MDS')
legend1 = ax[0].legend(*scatter.legend_elements(), title="Species")
ax[0].add_artist(legend1)

# NMDS plot
scatter = ax[1].scatter(X_nmds[:, 0], X_nmds[:, 1], c=y.astype('category').cat.codes)
ax[1].set_title('NMDS')
legend1 = ax[1].legend(*scatter.legend_elements(), title="Species")
ax[1].add_artist(legend1)

plt.show()
```

![output_3_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/211162/262a5fb2-f698-0267-da37-549476d8fe7c.png)


このコードでは、まずペンギンのデータセットをロードし、数値データだけを取り出しています（種類などのカテゴリデータは除く）。その後、MDSとNMDSをそれぞれのデータに適用し、結果をプロットしています。プロットの色はペンギンの種類によって異なります（`species`の値で色分け）。それぞれのプロットで、種類の違うペンギンがどのように分離するかをチェックしてみてください。

---

### **user**:

load_digitsのデータセットに対して、MDSとNMDSを実行し、結果を可視化して比較できるようなコードを書いてください。

---

### **assistant**:

Pythonのscikit-learnライブラリから手書き数字のデータセット`load_digits`を読み込み、そのデータに対してまずMDSを適用し、次にNMDSを適用して、その結果を比較するためのコードを以下に示します。

```python
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.manifold import MDS
from sklearn.metrics import pairwise_distances

# Load digits dataset
digits = datasets.load_digits()
X = digits.data
y = digits.target

# Apply MDS
mds = MDS(n_components=2, random_state=1)
X_mds = mds.fit_transform(X)

# Compute pairwise distance matrix
D = pairwise_distances(X)

# Apply non-metric MDS
nmds = MDS(n_components=2, dissimilarity='precomputed', random_state=1)
X_nmds = nmds.fit_transform(D)

# Plot results
fig, ax = plt.subplots(1, 2, figsize=(12, 6))

# MDS plot
ax[0].scatter(X_mds[:, 0], X_mds[:, 1], c=y)
ax[0].set_title('MDS')

# NMDS plot
ax[1].scatter(X_nmds[:, 0], X_nmds[:, 1], c=y)
ax[1].set_title('NMDS')

plt.show()
```

![output_4_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/211162/14d80de7-8d34-5312-c640-22d1e5057094.png)


このコードでは、まず手書き数字のデータセットをロードします。そのデータセットに対してMDSを適用し、次にNMDSを適用します。そして、その結果をプロットして可視化します。

プロットの結果は、手書き数字（0から9までの10個のクラスに対して、yの値で色分け）が2次元空間にどのように配置されているかを示しています。

---

### **user**:

load_linnerudのデータセットに対して、MDSとNMDSを実行し、結果を可視化して比較できるようなコードを書いてください。

---

### **assistant**:

Pythonのscikit-learnライブラリからリネルドのデータセット`load_linnerud`を読み込み、そのデータに対してまずMDSを適用し、次にNMDSを適用して、その結果を比較するためのコードを以下に示します。

なお、`load_linnerud`データセットは、フィットネスクラブの20人の中年男性のフィジカルフィットネステストに関するデータで、測定されたパラメータは３つ（体重、ウエスト、脈拍）です。

```python
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.manifold import MDS
from sklearn.metrics import pairwise_distances


# Load linnerud dataset
linnerud = datasets.load_linnerud()
X = linnerud.data
y = linnerud.target[:, 0]  #  Use the first target (Chin-ups) as labels for plotting

# Apply MDS
mds = MDS(n_components=2, random_state=1)
X_mds = mds.fit_transform(X)

# Compute pairwise distance matrix
D = pairwise_distances(X)

# Apply non-metric MDS
nmds = MDS(n_components=2, dissimilarity='precomputed', random_state=1)
X_nmds = nmds.fit_transform(D)

# Plot results
fig, ax = plt.subplots(1, 2, figsize=(12, 6))

# MDS plot
scatter = ax[0].scatter(X_mds[:, 0], X_mds[:, 1], c=y)
ax[0].set_title('MDS')
plt.colorbar(scatter, ax=ax[0])

# NMDS plot
scatter = ax[1].scatter(X_nmds[:, 0], X_nmds[:, 1], c=y)
ax[1].set_title('NMDS')
plt.colorbar(scatter, ax=ax[1])

plt.show()
```

![output_5_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/211162/d4ed622b-4a88-00df-67f4-454ab24a5977.png)


このコードでは、まずリネルドのデータセットをロードします。そのデータセットに対してMDSを適用し、次にNMDSを適用します。そして、その結果をプロットして可視化します。

プロットの色はチンアップ（懸垂）の回数（`y`の値）で、色の濃淡で示されています。これにより、各テスト参加者が2次元空間にどのように配置されているかを視覚的に確認することができます。

---

### **user**:

load_wineのデータセットに対して、MDSとNMDSを実行し、結果を可視化して比較できるようなコードを書いてください。

---

### **assistant**:

Pythonのscikit-learnライブラリからワインのデータセット`load_wine`を読み込み、そのデータに対してまずMDSを適用し、次にNMDSを適用して、その結果を比較するためのコードを以下に示します。

```python
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.manifold import MDS
from sklearn.metrics import pairwise_distances


# Load wine dataset
wine = datasets.load_wine()
X = wine.data
y = wine.target

# Apply MDS
mds = MDS(n_components=2, random_state=1)
X_mds = mds.fit_transform(X)

# Compute pairwise distance matrix
D = pairwise_distances(X)

# Apply non-metric MDS
nmds = MDS(n_components=2, dissimilarity='precomputed', random_state=1)
X_nmds = nmds.fit_transform(D)

# Plot results
fig, ax = plt.subplots(1, 2, figsize=(12, 6))

# MDS plot
scatter = ax[0].scatter(X_mds[:, 0], X_mds[:, 1], c=y)
ax[0].set_title('MDS')
legend1 = ax[0].legend(*scatter.legend_elements(), title="Class")
ax[0].add_artist(legend1)

# NMDS plot
scatter = ax[1].scatter(X_nmds[:, 0], X_nmds[:, 1], c=y)
ax[1].set_title('NMDS')
legend1 = ax[1].legend(*scatter.legend_elements(), title="Class")
ax[1].add_artist(legend1)

plt.show()
```

![output_6_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/211162/d3ee5851-a85d-b369-7ceb-83bba578484a.png)


このコードでは、まずワインのデータセットをロードします。そのデータセットに対してMDSを適用し、次にNMDSを適用します。そして、その結果をプロットして可視化します。

プロットの色はワインのクラス（`y`の値）で色分けされ、色の違いでワインのクラスを判別できます。ここでは、各ワインが2次元空間にどのように配置されているかを視覚的に確認することができます。

---

### **user**:

load_breast_cancerのデータセットに対して、MDSとNMDSを実行し、結果を可視化して比較できるようなコードを書いてください。

---

### **assistant**:

Pythonのscikit-learnライブラリから乳がんのデータセット`load_breast_cancer`を読み込み、そのデータに対してまずMDSを適用し、次にNMDSを適用して、その結果を比較するためのコードを以下に示します。

```python
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.manifold import MDS
from sklearn.metrics import pairwise_distances


# Load breast cancer dataset
cancer = datasets.load_breast_cancer()
X = cancer.data
y = cancer.target

# Apply MDS
mds = MDS(n_components=2, random_state=1)
X_mds = mds.fit_transform(X)

# Compute pairwise distance matrix
D = pairwise_distances(X)

# Apply non-metric MDS
nmds = MDS(n_components=2, dissimilarity='precomputed', random_state=1)
X_nmds = nmds.fit_transform(D)

# Plot results
fig, ax = plt.subplots(1, 2, figsize=(12, 6))

# MDS plot
scatter = ax[0].scatter(X_mds[:, 0], X_mds[:, 1], c=y)
ax[0].set_title('MDS')
legend1 = ax[0].legend(*scatter.legend_elements(), title="Class")
ax[0].add_artist(legend1)

# NMDS plot
scatter = ax[1].scatter(X_nmds[:, 0], X_nmds[:, 1], c=y)
ax[1].set_title('NMDS')
legend1 = ax[1].legend(*scatter.legend_elements(), title="Class")
ax[1].add_artist(legend1)

plt.show()
```

![output_7_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/211162/a78152bd-c5e4-ea2d-401b-ed12c6eb9105.png)


このコードでは、まず乳がんのデータセットをロードします。そのデータセットに対してMDSを適用し、次にNMDSを適用します。そして、その結果をプロットして可視化します。

プロットの色は乳がんのクラス（良性または悪性、`y`の値）で色分けされており、色の違いでクラスを判別できます。ここでは、各サンプルが2次元空間にどのように配置されているかを視覚的に確認することができます。

---

### **user**:

うーん、どれを見てもMDSとNMDSで違いがわからないですね。

---

### **assistant**:

確かに視覚的なプロットだけを見ていると、MDSとNMDSの違いは一見して分かりづらいことがあります。しかし、これら二つの手法は、データの距離関係をどのように扱うかにおいて基本的な違いがあります。

MDS（多次元尺度法）は、元の高次元データの距離関係をできるだけ忠実に低次元空間に再現しようと試みます。つまり、もとの高次元空間での距離が近いデータ点は、低次元空間でも距離が近くなるようにデータ点を配置します。

一方、NMDS（非尺度多次元尺度法）は、元の高次元データの距離の順序関係を保つことに重きを置きます。つまり、高次元空間で距離が近いデータのペアが、低次元空間でも距離が近いペアであることを保証しますが、その具体的な距離は必ずしも一致しないということです。

したがって、データにどの手法を適用するべきかを決定する際には、データの特性や、分析の目的を考慮することが重要です。また、結果の解釈においても、これらの手法の違いを理解しておくことが大切です。

---

### **user**:

pairwise distance matrix を違う種類の距離に変えたら変わりますかね？

---

### **assistant**:

はい、絶対に変わります。実際、MDSやNMDSで使用される距離計算は非常に重要な役割を果たします。

MDSやNMDSは距離行列を入力として使用するので、データ間の距離の計算方法が変わると、この距離行列自体が変わります。これにより、元の高次元データを低次元空間にマッピングする結果が異なることになります。

通常、デフォルトではユークリッド距離が使用されますが、他の距離メトリクスを使用することも可能です。たとえば、マンハッタン距離、マハラノビス距離、コサイン距離などがあります。これらの距離メトリクスは、データの特性や分析の目的によって選択することができます。

したがって、距離メトリクスを変えてみることで、データを異なる視点から解釈することが可能になります。また、これにより、異なるパターンや構造を見つけることができるかもしれません。

---

### **user**:

さまざまな距離メトリクスを変えてNMDSした結果を可視化して比較するコードを書いてください。

---

### **assistant**:

Pythonのsklearnおよびscipyライブラリを用いて、ユークリッド距離、マンハッタン距離、マハラノビス距離を計算し、それぞれの距離に基づいたNMDSの結果を比較するためのコードを以下に示します。

```python
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.manifold import MDS
from sklearn.metrics import pairwise_distances
from scipy.spatial.distance import pdist, squareform

# Load breast cancer dataset
cancer = datasets.load_breast_cancer()
X = cancer.data
y = cancer.target

# Compute pairwise distances
D_euclidean = pairwise_distances(X, metric='euclidean')
D_manhattan = pairwise_distances(X, metric='manhattan')
D_mahalanobis = squareform(pdist(X, metric='mahalanobis'))

# Apply NMDS with different distance metric
mds = MDS(n_components=2, dissimilarity='precomputed', random_state=1)
X_nmds_euclidean = mds.fit_transform(D_euclidean)
X_nmds_manhattan = mds.fit_transform(D_manhattan)
X_nmds_mahalanobis = mds.fit_transform(D_mahalanobis)

# Plot results
fig, ax = plt.subplots(1, 3, figsize=(18, 6))

# NMDS with Euclidean distance plot
scatter = ax[0].scatter(X_nmds_euclidean[:, 0], X_nmds_euclidean[:, 1], c=y)
ax[0].set_title('NMDS with Euclidean distance')
legend1 = ax[0].legend(*scatter.legend_elements(), title="Class")
ax[0].add_artist(legend1)

# NMDS with Manhattan distance plot
scatter = ax[1].scatter(X_nmds_manhattan[:, 0], X_nmds_manhattan[:, 1], c=y)
ax[1].set_title('NMDS with Manhattan distance')
legend1 = ax[1].legend(*scatter.legend_elements(), title="Class")
ax[1].add_artist(legend1)

# NMDS with Mahalanobis distance plot
scatter = ax[2].scatter(X_nmds_mahalanobis[:, 0], X_nmds_mahalanobis[:, 1], c=y)
ax[2].set_title('NMDS with Mahalanobis distance')
legend1 = ax[2].legend(*scatter.legend_elements(), title="Class")
ax[2].add_artist(legend1)

plt.show()
```

![output_8_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/211162/92a7af97-06cc-4401-6f48-ff5e2466e780.png)


このコードでは、ユークリッド距離、マンハッタン距離、マハラノビス距離を用いて距離行列を計算し、それぞれの結果を用いてNMDSを適用しています。その結果をそれぞれプロットして視覚化し、それらの違いを比較しています。

---

### **user**:

load_wineデータに対して、さまざまな距離メトリクスを変えてNMDSした結果を可視化して比較するコードを書いてください。

---

### **assistant**:

Pythonのsklearnおよびscipyライブラリを用いて、ユークリッド距離、マンハッタン距離、マハラノビス距離を計算し、それぞれの距離に基づいたNMDSの結果を比較するためのコードを以下に示します。ただし、マハラノビス距離を計算する際、逆共分散行列が必要となりますが、load_wineデータセットに対しては逆共分散行列を一意に求めることができないため、ここでは除外します。

```python
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.manifold import MDS
from sklearn.metrics import pairwise_distances

# Load wine dataset
wine = datasets.load_wine()
X = wine.data
y = wine.target

# Compute pairwise distances
D_euclidean = pairwise_distances(X, metric='euclidean')
D_manhattan = pairwise_distances(X, metric='manhattan')

# Apply NMDS with different distance metric
mds = MDS(n_components=2, dissimilarity='precomputed', random_state=1)
X_nmds_euclidean = mds.fit_transform(D_euclidean)
X_nmds_manhattan = mds.fit_transform(D_manhattan)

# Plot results
fig, ax = plt.subplots(1, 2, figsize=(12, 6))

# NMDS with Euclidean distance plot
scatter = ax[0].scatter(X_nmds_euclidean[:, 0], X_nmds_euclidean[:, 1], c=y)
ax[0].set_title('NMDS with Euclidean distance')
legend1 = ax[0].legend(*scatter.legend_elements(), title="Class")
ax[0].add_artist(legend1)

# NMDS with Manhattan distance plot
scatter = ax[1].scatter(X_nmds_manhattan[:, 0], X_nmds_manhattan[:, 1], c=y)
ax[1].set_title('NMDS with Manhattan distance')
legend1 = ax[1].legend(*scatter.legend_elements(), title="Class")
ax[1].add_artist(legend1)

plt.show()
```

![output_9_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/211162/9283cbf6-de0e-8538-ef14-240d1a04c0e7.png)


このコードでは、ユークリッド距離、マンハッタン距離を用いて距離行列を計算し、それぞれの結果を用いてNMDSを適用しています。その結果をそれぞれプロットして視覚化し、それらの違いを比較しています。

---

### **user**:

さすがに距離のメトリクスを変えたら変わりましたね

---

### **assistant**:

はい、距離メトリクスを変更すると、各データポイント間の距離関係が変わるため、その結果として可視化の結果も変わります。 

これは、選択した距離メトリクスがデータの特性をどのように捉えるかによって、元の高次元データを2次元へマッピングする結果が異なるためです。 

したがって、データの特性や解析の目的により、適切な距離メトリクスを選択することが重要となります。また、異なる距離メトリクスを用いて解析を行うことで、データの異なる側面を見ることができます。

---

### **user**:

NMDSについてより深く学ぶために、sklearn.manifold の MDS を使わずにNMDSをスクラッチ実装し、先ほどのデータをNMDSで可視化してください。

---

### **assistant**:

以下は、Scikit-Learnの`MDS`を使わないで、スクラッチでNMDSを行うための実装およびプロットの例です。ここでは、ワインのデータセット`load_wine`に対してスクラッチでNMDSを適用し、結果をプロットしています。

```python
import numpy as np
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.manifold import MDS
from sklearn.metrics import pairwise_distances
from scipy.optimize import minimize

def stress(X, D):
    D_hat = pairwise_distances(X)
    return ((D - D_hat)**2).sum()

# Load wine dataset
wine = datasets.load_wine()
X = wine.data
y = wine.target

n_samples = X.shape[0]
n_components = 2

# Compute pairwise distances
D = pairwise_distances(X, metric='euclidean')

# Initialize X randomly
X_init = np.random.RandomState(seed=0).randn(n_samples, n_components) 

# Minimize the stress
res = minimize(stress, X_init.flatten(), args=(D,), method='CG')

# The result (reshaped to the original form)
X_nmds_scratch = res.x.reshape(n_samples, n_components)

# Plot results
plt.scatter(X_nmds_scratch[:, 0], X_nmds_scratch[:, 1], c=y)
plt.title('NMDS with Euclidean distance (scratch)')
plt.gca().legend(*plt.gca().get_legend_handles_labels(), title="Class")
plt.show()
```

このスクラッチの実装では、まずワインのデータセットをロードし、そのデータに対してユークリッド距離を計算しています。そして、2次元の初期点をランダムに生成して、その中から`scipy.optimize.minimize`を使ってストレス関数を最小化する点を見つけています。

ストレス関数はMDSやNMDSで使用され、距離行列の元の距離`D`と低次元空間での距離`D_hat`の差の二乗和を表します。これを最小化することで、高次元空間の距離関係をできるだけ保ったまま低次元空間にデータをマッピングします。

最後に見つかった点をプロットして可視化します。このプロットでは、各ワインが2次元空間にどのように配置されているかを視覚的に確認することができます。

---

### **user**:

次のエラーを修正した全コードを書いてください。

---

### **assistant**:

すみませんが、具体的にどの部分がエラーを引き起こしているのかを教えていただくことは可能でしょうか？エラーとなる具体的なコードやエラーメッセージを提供いただけますと、より具体的なサポートが提供できます。

---

### **user**:

```
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[11], line 27
     24 X_init = np.random.RandomState(seed=0).randn(n_samples, n_components) 
     26 # Minimize the stress
---> 27 res = minimize(stress, X_init.flatten(), args=(D,), method='CG')
     29 # The result (reshaped to the original form)
     30 X_nmds_scratch = res.x.reshape(n_samples, n_components)

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_minimize.py:703, in minimize(fun, x0, args, method, jac, hess, hessp, bounds, constraints, tol, callback, options)
    701     res = _minimize_powell(fun, x0, args, callback, bounds, **options)
    702 elif meth == 'cg':
--> 703     res = _minimize_cg(fun, x0, args, jac, callback, **options)
    704 elif meth == 'bfgs':
    705     res = _minimize_bfgs(fun, x0, args, jac, callback, **options)

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_optimize.py:1753, in _minimize_cg(fun, x0, args, jac, callback, gtol, norm, eps, maxiter, disp, return_all, finite_diff_rel_step, **unknown_options)
   1750 if maxiter is None:
   1751     maxiter = len(x0) * 200
-> 1753 sf = _prepare_scalar_function(fun, x0, jac=jac, args=args, epsilon=eps,
   1754                               finite_diff_rel_step=finite_diff_rel_step)
   1756 f = sf.fun
   1757 myfprime = sf.grad

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_optimize.py:383, in _prepare_scalar_function(fun, x0, jac, args, bounds, epsilon, finite_diff_rel_step, hess)
    379     bounds = (-np.inf, np.inf)
    381 # ScalarFunction caches. Reuse of fun(x) during grad
    382 # calculation reduces overall function evaluations.
--> 383 sf = ScalarFunction(fun, x0, args, grad, hess,
    384                     finite_diff_rel_step, bounds, epsilon=epsilon)
    386 return sf

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py:158, in ScalarFunction.__init__(self, fun, x0, args, grad, hess, finite_diff_rel_step, finite_diff_bounds, epsilon)
    155     self.f = fun_wrapped(self.x)
    157 self._update_fun_impl = update_fun
--> 158 self._update_fun()
    160 # Gradient evaluation
    161 if callable(grad):

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py:251, in ScalarFunction._update_fun(self)
    249 def _update_fun(self):
    250     if not self.f_updated:
--> 251         self._update_fun_impl()
    252         self.f_updated = True

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py:155, in ScalarFunction.__init__.<locals>.update_fun()
    154 def update_fun():
--> 155     self.f = fun_wrapped(self.x)

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py:137, in ScalarFunction.__init__.<locals>.fun_wrapped(x)
    133 self.nfev += 1
    134 # Send a copy because the user may overwrite it.
    135 # Overwriting results in undefined behaviour because
    136 # fun(self.x) will change self.x, with the two no longer linked.
--> 137 fx = fun(np.copy(x), *args)
    138 # Make sure the function returns a true scalar
    139 if not np.isscalar(fx):

Cell In[11], line 9, in stress(X, D)
      8 def stress(X, D):
----> 9     D_hat = pairwise_distances(X)
     10     return ((D - D_hat)**2).sum()

File /usr/local/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:2196, in pairwise_distances(X, Y, metric, n_jobs, force_all_finite, **kwds)
   2193         return distance.squareform(distance.pdist(X, metric=metric, **kwds))
   2194     func = partial(distance.cdist, metric=metric, **kwds)
-> 2196 return _parallel_pairwise(X, Y, func, n_jobs, **kwds)

File /usr/local/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:1766, in _parallel_pairwise(X, Y, func, n_jobs, **kwds)
   1763 X, Y, dtype = _return_float_dtype(X, Y)
   1765 if effective_n_jobs(n_jobs) == 1:
-> 1766     return func(X, Y, **kwds)
   1768 # enforce a threading backend to prevent data communication overhead
   1769 fd = delayed(_dist_wrapper)

File /usr/local/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:310, in euclidean_distances(X, Y, Y_norm_squared, squared, X_norm_squared)
    234 def euclidean_distances(
    235     X, Y=None, *, Y_norm_squared=None, squared=False, X_norm_squared=None
    236 ):
    237     """
    238     Compute the distance matrix between each pair from a vector array X and Y.
    239 
   (...)
    308            [1.41421356]])
    309     """
--> 310     X, Y = check_pairwise_arrays(X, Y)
    312     if X_norm_squared is not None:
    313         X_norm_squared = check_array(X_norm_squared, ensure_2d=False)

File /usr/local/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:156, in check_pairwise_arrays(X, Y, precomputed, dtype, accept_sparse, force_all_finite, copy)
    153     dtype = dtype_float
    155 if Y is X or Y is None:
--> 156     X = Y = check_array(
    157         X,
    158         accept_sparse=accept_sparse,
    159         dtype=dtype,
    160         copy=copy,
    161         force_all_finite=force_all_finite,
    162         estimator=estimator,
    163     )
    164 else:
    165     X = check_array(
    166         X,
    167         accept_sparse=accept_sparse,
   (...)
    171         estimator=estimator,
    172     )

File /usr/local/lib/python3.10/site-packages/sklearn/utils/validation.py:938, in check_array(array, accept_sparse, accept_large_sparse, dtype, order, copy, force_all_finite, ensure_2d, allow_nd, ensure_min_samples, ensure_min_features, estimator, input_name)
    936     # If input is 1D raise error
    937     if array.ndim == 1:
--> 938         raise ValueError(
    939             "Expected 2D array, got 1D array instead:\narray={}.\n"
    940             "Reshape your data either using array.reshape(-1, 1) if "
    941             "your data has a single feature or array.reshape(1, -1) "
    942             "if it contains a single sample.".format(array)
    943         )
    945 if dtype_numeric and hasattr(array.dtype, "kind") and array.dtype.kind in "USV":
    946     raise ValueError(
    947         "dtype='numeric' is not compatible with arrays of bytes/strings."
    948         "Convert your data to numeric values explicitly instead."
    949     )

ValueError: Expected 2D array, got 1D array instead:
array=[ 1.76405235  0.40015721  0.97873798  2.2408932   1.86755799 -0.97727788
  0.95008842 -0.15135721 -0.10321885  0.4105985   0.14404357  1.45427351
  0.76103773  0.12167502  0.44386323  0.33367433  1.49407907 -0.20515826
  0.3130677  -0.85409574 -2.55298982  0.6536186   0.8644362  -0.74216502
  2.26975462 -1.45436567  0.04575852 -0.18718385  1.53277921  1.46935877
  0.15494743  0.37816252 -0.88778575 -1.98079647 -0.34791215  0.15634897
  1.23029068  1.20237985 -0.38732682 -0.30230275 -1.04855297 -1.42001794
 -1.70627019  1.9507754  -0.50965218 -0.4380743  -1.25279536  0.77749036
 -1.61389785 -0.21274028 -0.89546656  0.3869025  -0.51080514 -1.18063218
 -0.02818223  0.42833187  0.06651722  0.3024719  -0.63432209 -0.36274117
 -0.67246045 -0.35955316 -0.81314628 -1.7262826   0.17742614 -0.40178094
 -1.63019835  0.46278226 -0.90729836  0.0519454   0.72909056  0.12898291
  1.13940068 -1.23482582  0.40234164 -0.68481009 -0.87079715 -0.57884966
 -0.31155253  0.05616534 -1.16514984  0.90082649  0.46566244 -1.53624369
  1.48825219  1.89588918  1.17877957 -0.17992484 -1.07075262  1.05445173
 -0.40317695  1.22244507  0.20827498  0.97663904  0.3563664   0.70657317
  0.01050002  1.78587049  0.12691209  0.40198936  1.8831507  -1.34775906
 -1.270485    0.96939671 -1.17312341  1.94362119 -0.41361898 -0.74745481
  1.92294203  1.48051479  1.86755896  0.90604466 -0.86122569  1.91006495
 -0.26800337  0.8024564   0.94725197 -0.15501009  0.61407937  0.92220667
  0.37642553 -1.09940079  0.29823817  1.3263859  -0.69456786 -0.14963454
 -0.43515355  1.84926373  0.67229476  0.40746184 -0.76991607  0.53924919
 -0.67433266  0.03183056 -0.63584608  0.67643329  0.57659082 -0.20829876
  0.39600671 -1.09306151 -1.49125759  0.4393917   0.1666735   0.63503144
  2.38314477  0.94447949 -0.91282223  1.11701629 -1.31590741 -0.4615846
 -0.06824161  1.71334272 -0.74475482 -0.82643854 -0.09845252 -0.66347829
  1.12663592 -1.07993151 -1.14746865 -0.43782004 -0.49803245  1.92953205
  0.94942081  0.08755124 -1.22543552  0.84436298 -1.00021535 -1.5447711
  1.18802979  0.31694261  0.92085882  0.31872765  0.85683061 -0.65102559
 -1.03424284  0.68159452 -0.80340966 -0.68954978 -0.4555325   0.01747916
 -0.35399391 -1.37495129 -0.6436184  -2.22340315  0.62523145 -1.60205766
 -1.10438334  0.05216508 -0.739563    1.5430146  -1.29285691  0.26705087
 -0.03928282 -1.1680935   0.52327666 -0.17154633  0.77179055  0.82350415
  2.16323595  1.33652795 -0.36918184 -0.23937918  1.0996596   0.65526373
  0.64013153 -1.61695604 -0.02432612 -0.73803091  0.2799246  -0.09815039
  0.91017891  0.31721822  0.78632796 -0.4664191  -0.94444626 -0.41004969
 -0.01702041  0.37915174  2.25930895 -0.04225715 -0.955945   -0.34598178
 -0.46359597  0.48148147 -1.54079701  0.06326199  0.15650654  0.23218104
 -0.59731607 -0.23792173 -1.42406091 -0.49331988 -0.54286148  0.41605005
 -1.15618243  0.7811981   1.49448454 -2.06998503  0.42625873  0.67690804
 -0.63743703 -0.39727181 -0.13288058 -0.29779088 -0.30901297 -1.67600381
  1.15233156  1.07961859 -0.81336426 -1.46642433  0.52106488 -0.57578797
  0.14195316 -0.31932842  0.69153875  0.69474914 -0.72559738 -1.38336396
 -1.5829384   0.61037938 -1.18885926 -0.50681635 -0.59631404 -0.0525673
 -1.93627981  0.1887786   0.52389102  0.08842209 -0.31088617  0.09740017
  0.39904635 -2.77259276  1.95591231  0.39009332 -0.65240858 -0.39095338
  0.49374178 -0.11610394 -2.03068447  2.06449286 -0.11054066  1.02017271
 -0.69204985  1.53637705  0.28634369  0.60884383 -1.04525337  1.21114529
  0.68981816  1.30184623 -0.62808756 -0.48102712  2.3039167  -1.06001582
 -0.1359497   1.13689136  0.09772497  0.58295368 -0.39944903  0.37005589
 -1.30652685  1.65813068 -0.11816405 -0.6801782   0.66638308 -0.46071979
 -1.33425847 -1.34671751  0.69377315 -0.15957344 -0.13370156  1.07774381
 -1.12682581 -0.73067775 -0.38487981  0.09435159 -0.04217145 -0.28688719
 -0.0616264  -0.10730528 -0.71960439 -0.81299299  0.27451636 -0.89091508
 -1.15735526 -0.31229225 -0.15766702  2.2567235  -0.70470028  0.94326072
  0.74718833 -1.18894496  0.77325298 -1.18388064 -2.65917224  0.60631952
 -1.75589058  0.45093446 -0.6840109   1.6595508   1.0685094  -0.4533858
 -0.68783761 -1.2140774  -0.44092263 -0.2803555  -0.36469354  0.15670386
  0.5785215   0.34965446 -0.76414392 -1.43779147  1.36453185 -0.68944918
 -0.6522936  -0.52118931].
Reshape your data either using array.reshape(-1, 1) if your data has a single feature or array.reshape(1, -1) if it contains a single sample.
```

---

### **assistant**:



---

### **user**:

```
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[11], line 27
     24 X_init = np.random.RandomState(seed=0).randn(n_samples, n_components) 
     26 # Minimize the stress
---> 27 res = minimize(stress, X_init.flatten(), args=(D,), method='CG')
     29 # The result (reshaped to the original form)
     30 X_nmds_scratch = res.x.reshape(n_samples, n_components)

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_minimize.py:703, in minimize(fun, x0, args, method, jac, hess, hessp, bounds, constraints, tol, callback, options)
    701     res = _minimize_powell(fun, x0, args, callback, bounds, **options)
    702 elif meth == 'cg':
--> 703     res = _minimize_cg(fun, x0, args, jac, callback, **options)
    704 elif meth == 'bfgs':
    705     res = _minimize_bfgs(fun, x0, args, jac, callback, **options)

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_optimize.py:1753, in _minimize_cg(fun, x0, args, jac, callback, gtol, norm, eps, maxiter, disp, return_all, finite_diff_rel_step, **unknown_options)
   1750 if maxiter is None:
   1751     maxiter = len(x0) * 200
-> 1753 sf = _prepare_scalar_function(fun, x0, jac=jac, args=args, epsilon=eps,
   1754                               finite_diff_rel_step=finite_diff_rel_step)
   1756 f = sf.fun
   1757 myfprime = sf.grad

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_optimize.py:383, in _prepare_scalar_function(fun, x0, jac, args, bounds, epsilon, finite_diff_rel_step, hess)
    379     bounds = (-np.inf, np.inf)
    381 # ScalarFunction caches. Reuse of fun(x) during grad
    382 # calculation reduces overall function evaluations.
--> 383 sf = ScalarFunction(fun, x0, args, grad, hess,
    384                     finite_diff_rel_step, bounds, epsilon=epsilon)
    386 return sf

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py:158, in ScalarFunction.__init__(self, fun, x0, args, grad, hess, finite_diff_rel_step, finite_diff_bounds, epsilon)
    155     self.f = fun_wrapped(self.x)
    157 self._update_fun_impl = update_fun
--> 158 self._update_fun()
    160 # Gradient evaluation
    161 if callable(grad):

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py:251, in ScalarFunction._update_fun(self)
    249 def _update_fun(self):
    250     if not self.f_updated:
--> 251         self._update_fun_impl()
    252         self.f_updated = True

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py:155, in ScalarFunction.__init__.<locals>.update_fun()
    154 def update_fun():
--> 155     self.f = fun_wrapped(self.x)

File /usr/local/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py:137, in ScalarFunction.__init__.<locals>.fun_wrapped(x)
    133 self.nfev += 1
    134 # Send a copy because the user may overwrite it.
    135 # Overwriting results in undefined behaviour because
    136 # fun(self.x) will change self.x, with the two no longer linked.
--> 137 fx = fun(np.copy(x), *args)
    138 # Make sure the function returns a true scalar
    139 if not np.isscalar(fx):

Cell In[11], line 9, in stress(X, D)
      8 def stress(X, D):
----> 9     D_hat = pairwise_distances(X)
     10     return ((D - D_hat)**2).sum()

File /usr/local/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:2196, in pairwise_distances(X, Y, metric, n_jobs, force_all_finite, **kwds)
   2193         return distance.squareform(distance.pdist(X, metric=metric, **kwds))
   2194     func = partial(distance.cdist, metric=metric, **kwds)
-> 2196 return _parallel_pairwise(X, Y, func, n_jobs, **kwds)

File /usr/local/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:1766, in _parallel_pairwise(X, Y, func, n_jobs, **kwds)
   1763 X, Y, dtype = _return_float_dtype(X, Y)
   1765 if effective_n_jobs(n_jobs) == 1:
-> 1766     return func(X, Y, **kwds)
   1768 # enforce a threading backend to prevent data communication overhead
   1769 fd = delayed(_dist_wrapper)

File /usr/local/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:310, in euclidean_distances(X, Y, Y_norm_squared, squared, X_norm_squared)
    234 def euclidean_distances(
    235     X, Y=None, *, Y_norm_squared=None, squared=False, X_norm_squared=None
    236 ):
    237     """
    238     Compute the distance matrix between each pair from a vector array X and Y.
    239 
   (...)
    308            [1.41421356]])
    309     """
--> 310     X, Y = check_pairwise_arrays(X, Y)
    312     if X_norm_squared is not None:
    313         X_norm_squared = check_array(X_norm_squared, ensure_2d=False)

File /usr/local/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:156, in check_pairwise_arrays(X, Y, precomputed, dtype, accept_sparse, force_all_finite, copy)
    153     dtype = dtype_float
    155 if Y is X or Y is None:
--> 156     X = Y = check_array(
    157         X,
    158         accept_sparse=accept_sparse,
    159         dtype=dtype,
    160         copy=copy,
    161         force_all_finite=force_all_finite,
    162         estimator=estimator,
    163     )
    164 else:
    165     X = check_array(
    166         X,
    167         accept_sparse=accept_sparse,
   (...)
    171         estimator=estimator,
    172     )

File /usr/local/lib/python3.10/site-packages/sklearn/utils/validation.py:938, in check_array(array, accept_sparse, accept_large_sparse, dtype, order, copy, force_all_finite, ensure_2d, allow_nd, ensure_min_samples, ensure_min_features, estimator, input_name)
    936     # If input is 1D raise error
    937     if array.ndim == 1:
--> 938         raise ValueError(
    939             "Expected 2D array, got 1D array instead:\narray={}.\n"
    940             "Reshape your data either using array.reshape(-1, 1) if "
    941             "your data has a single feature or array.reshape(1, -1) "
    942             "if it contains a single sample.".format(array)
    943         )
    945 if dtype_numeric and hasattr(array.dtype, "kind") and array.dtype.kind in "USV":
    946     raise ValueError(
    947         "dtype='numeric' is not compatible with arrays of bytes/strings."
    948         "Convert your data to numeric values explicitly instead."
    949     )

ValueError: Expected 2D array, got 1D array instead:
Reshape your data either using array.reshape(-1, 1) if your data has a single feature or array.reshape(1, -1) if it contains a single sample.
```

---

### **assistant**:

This error is occurring because the function `pairwise_distances` from `sklearn.metrics` is expecting a 2D array but is receiving a 1D array. The error message suggests to reshape your data either using `array.reshape(-1, 1)` if your data has a single feature or `array.reshape(1, -1)` if it contains a single sample.

On line 27, when you call `minimize(stress, X_init.flatten(), args=(D,), method='CG')`, `X_init.flatten()` flattens your `X_init` data into 1D. Then, inside your `stress` function, you're passing this 1D data into `pairwise_distances`, which is throwing the error. 

You should be able to fix this error by reshaping the data back to 2D before passing it to 'pairwise_distances'. Here's a way to do that:

```python
def stress(X, D):
    X = X.reshape(-1, 2)  # Reshape back to 2D
    D_hat = pairwise_distances(X)
    return ((D - D_hat)**2).sum()
```

This will reshape your data back into a 2D array with 2 columns (the number of features in your original data) before calling `pairwise_distances`.

![output_12_0.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/211162/9dc592f4-907b-1324-be0c-f2edfe1b5a19.png)


---

### **user**:

できました。でもけっこう時間かかりましたね。

---

### **assistant**:

はい、最適化の計算は特に大きなデータセットでは時間がかかることがあります。計算時間を減らすためのいくつかの方法があります：

1. **初期値の選択**：最適化の初期値が適切でないと、アルゴリズムは最適な解に到達するのに時間がかかるか、全く到達しない可能性があります。異なる初期値を試すか、プロブレムに適した方法で初期値を選択すると良いでしょう。

2. **アルゴリズムの選択**：`minimize`関数は様々な最適化アルゴリズムを提供しています。`CG`メソッドは一般的には良い結果を提供しますが、あなたの特定の問題に対しては他のアルゴリズム（`BFGS`、`Newton-CG`、`L-BFGS-B`など）の方が効果的かもしれません。

3. 大規模な最適化問題のためのアルゴリズム（例えば、確率的勾配降下法）を使用すると、より速く解を見つけることができることがあります。

どの方法が最も効果的かは、具体的な問題と使用しているデータによります。一部の試行錯誤が必要となるかもしれません。

---

### **user**:

上記のスクラッチコードの中で、NMDSとMDSの違いを表す部分はどこでしょうか？

---

### **assistant**:

上記のコードでは、`Non-metric Multidimensional Scaling (NMDS)`と`Metric Multidimensional Scaling (MDS)`の違いは主に`stress`関数に表れます。具体的には以下の部分です。

```python
def stress(X, D):
    X = X.reshape(-1, 2)
    D_hat = pairwise_distances(X)
    return ((D - D_hat)**2).sum()
```
この関数は**Squared Stress**と呼ばれ、それは元の距離`D`（pairwise distance matrix）と再構成された距離`D_hat`(pairwise distance matrix of transformed data）の間の二乗誤差を合計したものです。

非メトリック多次元尺度分析（NMDS）では、このエネルギー関数を最小化することで、元のサンプル間の距離関係を尽可能維持した上でデータの次元を低減します。

一方、メトリック多次元尺度分析（MDS）では、`D`と`D_hat`の間の相関を最大化しようとします。非メトリックとメトリックのMDSの違いは、距離または類似性をどのようにモデル化するかによります。非メトリックMDSは、順序またはランクに基づく距離尺度を保持することに焦点を当て（つまり、距離の大小関係だけを考慮し、その絶対値は無視します）、メトリックMDSは実際の距離尺度を保持することに焦点を当てます。

ただし注意点として、上記のコードでは厳密なNMDSの実装ではありません。真のNMDSはデータ間の順序関係（距離のランク）を保存しようとしますが、このコードでは実際の距離を保持しようとしています。これはむしろメトリックMDSに近い形となります。

---

### **user**:

ええっと、すみません、結局、上記のコードはMDSなんですか？NMDSなんですか？

---

### **assistant**:

申し訳ありません、説明が混乱させるものでした。

上記のコードは、特に距離の二乗誤差を最小化する目的関数（stress関数）を使用していますので、これは一般に`Metric Multidimensional Scaling (MDS)`に最も近いです。

`Non-metric Multidimensional Scaling (NMDS)`は、元の距離データの順序（ランク）だけを保持しようとするため、アルゴリズム自体と目的関数が異なります。具体的には、それぞれの距離ペア間の順序関係（一つが他のものよりも大きいまたは小さい）だけが結果に影響を与えます。

したがって、上記のコードは厳密には`MDS`に分類されます。しかし、これらのアルゴリズムの境界はかなり流動的であることを理解しておくことが重要です。`MDS`および`NMDS`の両方がデータの次元削減を試み、元のデータ空間の構造と距離を保つことを試みますが、その方法はアルゴリズムの詳細と目的関数の選択により異なります。

---

### **user**:

上記のstrees関数を、MDSではなくHMDSになるように改変したコードを書いてください。

---

### **assistant**:

標準的な `Non-metric Multidimensional Scaling (NMDS)` アプローチは Kruskal の Rank-1 Stress を最小化することで、データポイント間の距離のランクを保存します。ただし、この目的関数は非常に複雑で、通常は多次元尺度法の実装に組み込まれた最適化ルーチンを使用して解決されます。

このため、一般的な最適化関数 `minimize` でこの問題を解決するのは難しいかもしれません。Python向けのNMDSアルゴリズムを実装するためのライブラリは、 `sklearn.manifold.MDS` で提供されており、ここでは `n_components` と `metric` というパラメータを使用します。

以下に具体的なNMDSを実行するサンプルコードを記載します：

```python
from sklearn.manifold import MDS
from sklearn.metrics import pairwise_distances

# Compute Pairwise distances for your matrix 'data'
D = pairwise_distances(data)

# Non-metric MDS
nmds = MDS(n_components=2, metric=False, dissimilarity='precomputed')
X_nmds_nonmetric = nmds.fit_transform(D)
```

このコードでは、`metric=False` とすることでNMDSを有効にし、`dissimilarity='precomputed'` として既に計算済みの距離行列を利用し、距離行列のランクの保存を試みます。

しかし、もしstrees関数を実装したいということであれば、以下のような関数を考えることができます。ただし、これはあくまで参考の一つであり、正確な形状は問題や扱っている具体的なデータによって変わります。

```python
def rank_order_stress(X, D):
    X = X.reshape(-1, 2)
    D_hat = pairwise_distances(X)
    
    # Convert distance matrices to rank matrices
    D = D.argsort().argsort()
    D_hat = D_hat.argsort().argsort()
    
    return ((D - D_hat)**2).sum()
```

![output_14_0.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/211162/c095403a-ef9f-b43b-cebf-167f124793f1.png)


しかし、このような非メトリックなストレス関数を用いて直接 'minimize' 関数を計算するのはあまり一般的ではありません。


---

### **user**:

計算に時間がかかるかと思ったら一瞬で答えが返ってきましたね。計算が正しいかどうかは微妙な気がしますが。

---

### **assistant**:

`scikit-learn`での多次元尺度法の実装は、特定の最適化ルーチンや初期化ストラテジなど、高速な計算を可能にする様々なテクニックを使用しています。そのため、計算が素早く行われることは一般的です。

しかしながら、結果が期待通りのものであるかどうかを確認するには、結果を可視化するか、別の方法（他の尺度法や主成分分析など）で同様の結果が得られるかどうかを確認すると良いでしょう。

多次元尺度法や他の次元削減法は探索的な分析手法であり、結果は常に適切な解釈と注意深い評価を必要とします。形状や配置がオリジナルデータの構造を適切に表現しているかどうかを評価するのに、領域の知識やデータの理解が重要となります。

---

### **user**:

まあいいや。だいたい分かりました。ありがとうございます。それではNMDSについてまとめの文章を詳しく書いてください。

---

### **assistant**:

承知しました。以下に`Non-metric Multidimensional Scaling（NMDS）`の詳細についての要約を記述いたします。

`Non-metric Multidimensional Scaling (NMDS)`は、高次元データの次元削減を試みる一連の技術の一つであり、データの構造やパターンを可視化するためによく用いられます。NMDSの目的は、高次元データ空間の順序（ランク）関係を保存することで、データの構造をなるべく損ねることなく二次元または三次元空間でデータを表現することです。

NMDSは、各データポイント間の距離の順序（ランク）だけを保存しようとします。距離の順序関係とは、あるデータポイントが他のデータポイントよりも近く（または遠く）にあるという事実です。この方法で、高次元空間でのデータの配置を低次元空間に「マッピング」することができます。

具体的な計算プロセスは以下の通りです：

1. 元の高次元データからペアワイズ距離の行列（または類似性行列）を計算します。
2. ランダムまたは特定の方法で低次元空間でのデータポイントの初期配置を生成します。
3. 低次元空間での各データポイント間のペアワイズ距離を計算します。
4. 低次元空間での距離と元の高次元空間での距離の順序を比較し、それらの間の「ストレス」を計算します。ストレスは、これら2つの距離配置がどれだけ一致していないかを示す指標であり、一般的には距離のランクの違いの平方和として計算されます。
5. ストレス関数を最小化するように低次元空間でのデータポイントの配置を最適化します。これは一般的に反復的な最適化手法で行われます。

Pythonでは、scikit-learnライブラリの`sklearn.manifold.MDS`クラスがNMDSを実装しています。ただし、厳密なNMDSの計算は計算的に複雑であり、結果の解釈は領域の知識とデータへの理解を必要とします。

---

### **user**:

最後に、以上の議論を踏まえて、MDSとNMDSに関する面白ギャグを創作してください。

---

### **assistant**:

それではいかがでしょう？

MDSとNMDSがバーで話しているところを想像してみてください。

MDSはNMDSに言います：「おい、私たちあんまり違わないよな。どちらも高次元のデータを低次元にプロジェクトするんだから。」

それに対して、NMDSは言いました：「たしかにそうだね。でも、ちょっと違うところがあるんだ。たとえば、君は狂ったように正確な距離を集計しようとするけど、僕はもっとリラックスして、許容範囲内で順序だけを維持しようとするんだよ。君が考えすぎだって言うことがあるけど、まあ、それが君のやり方だね。」

すると、MDSは言いました：「うーん、それなら次にデートに行くのは僕が招待するよ。君が今度は距離を気にしなくていいように。」

---