三点シーソーモデルについて (c.f., 伊藤 2004)

In [1]:
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
from mpl_toolkits.mplot3d import Axes3D
import random
In [2]:
m = 2000
#r1 = [random.normalvariate(mu=0,sigma=1.0) for i in range(0,n)]
#r2 = [random.normalvariate(mu=0,sigma=1.0) for i in range(0,n)]

factor = np.sqrt(3)*2.0 #標準偏差が1になるように規格化
r1 = 1.5*factor*np.array([random.random()-0.5  for i in range(0,m)]) 
r2 = factor*np.array([random.random()-0.5  for i in range(0,m)])
In [3]:
%matplotlib notebook

# (1) AOが真(Eq. 9.1)
xa = -r1 + r2 #x軸
xp = -r1 - r2 #y軸
xn = 2*r1 #z軸
xanom = np.vstack([xa,xp,xn]).T
print(xanom.shape)

# (2) AOが偽(NAO-PNA) (Eq. 9.2)
ya = -r1 
yp = -r2
yn = r1 + r2
yanom = np.vstack([ya,yp,yn]).T

# 特異値分解
xU, xs, xVh = np.linalg.svd(xanom)
yU, ys, yVh = np.linalg.svd(yanom)


# 描画
fig = plt.figure(figsize=(10,5))
ax = fig.add_subplot(121,projection='3d')
ax2 = fig.add_subplot(122,projection='3d')

ax.scatter(xa,xp,xn,alpha=0.3,s=5)
ax2.scatter(ya,yp,yn,alpha=0.3,s=5)

for i in range(0,2):
    x2 = xs[i]*xVh[i,0]/m**0.5; y2 = xs[i]*xVh[i,1]/m**0.5; z2 = xs[i]*xVh[i,2]/m**0.5
    ax.quiver([0],[0],[0],[x2],[y2],[z2],color='red') #主成分ベクトル(AO真)
    
    print("EOF" + str(i+1) +" for AO-NCM: ",[x2,y2,z2])
    
    x2 = ys[i]*yVh[i,0]/m**0.5; y2 = ys[i]*yVh[i,1]/m**0.5; z2 = ys[i]*yVh[i,2]/m**0.5
    ax2.quiver([0],[0],[0],[x2],[y2],[z2],color='red') #主成分ベクトル(AO偽)

    
    print("EOF" + str(i+1) +" for NAO-PNA: ",[x2,y2,z2])

    
# AO, PNAの基底ベクトル
ax.quiver([0],[0],[0],[-1],[-1],[2],color='orange')
ax.quiver([0],[0],[0],[1],[-1],[0],color='orange')
    
    
# NAO, PNAの基底ベクトル
ax2.quiver([0],[0],[0],[-1],[0],[1],color='orange')
ax2.quiver([0],[0],[0],[0],[-1],[1],color='orange')

ax.set_xlabel('Atlantic (A)'); ax.set_ylabel('Pacific (P)'); ax.set_zlabel('North Pole (N)')
ax2.set_xlabel('Atlantic (A)'); ax2.set_ylabel('Pacific (P)'); ax2.set_zlabel('North Pole (N)')

ax.set_xlim(-4,4);ax.set_ylim(-4,4);ax.set_zlim(-4,4)
ax2.set_xlim(-4,4);ax2.set_ylim(-4,4);ax2.set_zlim(-4,4)

plt.show()
(2000, 3)
EOF1 for AO-NCM:  [-1.479526442587139, -1.4968612601824351, 2.976387702769572]
EOF1 for NAO-PNA:  [-1.3453717205595774, -0.4330551386392751, 1.778426859198853]
EOF2 for AO-NCM:  [0.9994043003317551, -0.9955313901955901, -0.003872910136164721]
EOF2 for NAO-PNA:  [0.6361600605793124, -0.8985991727173677, 0.26243911213805493]

  • (左)AOが"真"の場合、主成分ベクトル(赤)は真のモードの基底ベクトル(オレンジ)と一致する。
  • (右)AOが"偽"の場合でも、主成分ベクトル(赤)はAOのようなパターンを示す(つまり、A-PとNで逆相関)。これは真のモードの基底ベクトル(オレンジ:NAO, PNA)と一致しない。

つまり、

三点データのEOF解析だけからでは何が真のモードなのか判断できない。

ことが分かる。


二点EOF(中緯度EOFに相当する)

In [4]:
xanom = np.vstack([xa,xp]).T

yanom = np.vstack([ya,yp]).T

# 特異値分解(Numpyの関数を使う)
xU, xs, xVh = np.linalg.svd(xanom)
yU, ys, yVh = np.linalg.svd(yanom)


# 描画
fig = plt.figure(figsize=(10,5))
ax = fig.add_subplot(121)
ax2 = fig.add_subplot(122)

ax.scatter(xa,xp,s=10)
ax2.scatter(ya,yp,s=10)

for i in range(0,2):
    x2 = xs[i]*xVh[i,0]/m**0.5; y2 = xs[i]*xVh[i,1]/m**0.5
    ax.quiver([0],[0],[x2],[y2],color='red',label ="EOF"+str(i))
    print("EOF" + str(i+1) +" for AO-NCM: ",[x2,y2])
    
    x2 = ys[i]*yVh[i,0]/m**0.5; y2 = ys[i]*yVh[i,1]/m**0.5
    ax2.quiver([0],[0],[x2],[y2],color='red',label ="EOF"+str(i))
    print("EOF" + str(i+1) +" for NAO-PNA: ",[x2,y2])
    
# NAO, PNAの基底ベクトル
ax.quiver([0],[0],[-1],[-1],color='orange')
ax.quiver([0],[0],[1],[-1],color='orange')
    
# NAO, PNAの基底ベクトル
ax2.quiver([0],[0],[-1],[0],color='orange')
ax2.quiver([0],[0],[0],[-1],color='orange')


ax.set_xlabel('Atlantic (A)'); ax.set_ylabel('Pacific (P)')
ax2.set_xlabel('Atlantic (A)'); ax2.set_ylabel('Pacific (P)')

ax.set_xlim(-4,4);ax.set_ylim(-4,4)
ax2.set_xlim(-4,4);ax2.set_ylim(-4,4)

plt.show()
EOF1 for AO-NCM:  [-1.4747880254614545, -1.5015481578252263]
EOF1 for NAO-PNA:  [-1.4881680916433402, -0.013380066181885841]
EOF2 for AO-NCM:  [-1.0063834904208544, 0.9884480314267304]
EOF2 for NAO-PNA:  [-0.008967729497061946, 0.9974157609237922]

 二点EOFの場合は、AOが"真"の場合(左)でも"偽"の場合(右)でも

主成分ベクトル(赤)が真のモードの基底ベクトル(オレンジ)と一致する

ことがわかる(符号の不定性は除いて)。

  • (左)AOが"真":主成分ベクトルは、PとAで相関(or 逆相関)を持つモード
  • (右)AOが"偽":主成分ベクトルは、PとAで独立なモード(つまり、NAOとPNAの基底ベクトル)。

 これが、AOの実在性を調べる為に中緯度EOFを施す根拠となっている。

In [ ]: