都道府県の長さ

ACTIVE GALATTICさんが，次のように書いていた。

中核都市の金沢と被害が激しい輪島は東京と水戸くらい離れていると石川県の広さを説明する投稿を先日見かけたが、石川県の広さは34~35位とランキング下位なので、むしろ東京や大阪が都道府県としては異様に狭いと認識すべきなのだろう

面積と距離が微妙に絡まっていて，ややモニョる。ある図形が与えられたときに，その面積ではなくて，長さを図形の大きさの指標として与えることはできるだろうか。地図上の領域だと東西南北の端点を考えて，南北間距離，東西間距離を抽出できるが，これでは2つの数字の組になる。あるいは図形の周長も考えられるが，フラクタル的な図形だと面積が等しくても周長はいくらでも長くできる。

そこで考えたのが次のような指標だ。図形を格子状の正方形に分割する。これらの正方形の全ての対を考え，各対の正方形の中心間の距離の平均値を求める。分割する正方形のサイズを小さくした極限として，一つの長さが定まるのではないか。（方法１）

検索してみると，もっとスマートな方法があった。図形内の点をランダムに選ぶ。その際に縦横方向にそれぞれ座標値が一様分布であるとした確率分布を仮定する。こうして選んだ2点間の距離の期待値を求めればよいというものだ。（方法２）

まず，方法１のプログラムをJuliaで作成して計算してみた。

a=zeros(Int64,1000001,2)

function ju(a,m,n)

  k = 0

  for i in 1:m

    for j in 1:n

      k = k + 1

      a[k,1] = i

      a[k,2] = j

    end

  end

end

function su(a,m,n)

  mn = m*n

  sum = 0

  for i = 1:mn

    for j = i:mn

      sum = sum + sqrt((a[i,1]-a[j,1])^2+(a[i,2]-a[j,2])^2)

    end

  end

  return sum/binomial(mn,2)

end

m=800

n=800

ju(a,m,n)

@time su(a,m,n)/m

n=   25    time= 0.016720　0.521829

n=   50　time= 0.018188　0.521513

n= 100　time= 0.065998　0.521433

n= 200　time= 0.788281　0.521412

n= 400　time= 12.1924  　0.521407

n= 800　time= 198.615  　0.521405

nは正方形の1辺の格子点の数，timeは計算時間(秒)であり，最後が計算結果の平均長である。

方法２による理論値が解析的に得られており，( 2+√2+5 log(√2 + 1) )/15 = 0.521405 なので，nを増やすとともに理論値に収束し，n=800では理論値に一致している。

m = 50; n = 50;

a = Flatten[Table[{i, j}, {i, 1, m}, {j, 1, n}], 1];

Timing[Sum[

     Sqrt[(a[[i, 1]] - a[[j, 1]])^2 + (a[[i, 2]] - a[[j, 2]])^2], 

　　　{i, 1, n*n}, {j, i, n*n}]/Binomial[n*n, 2]/n // N]

Mathematicaでは，若干コードが簡単になるけれど，時間は n=50 で 12.4801秒と，600倍以上かかってしまった。

（付）このコードはそのまま乱数による計算に置き換えることもできる。Juliaでは，a[k,1] = rand()，a[k,2] = rand() として最後の結果の/nを除く。Mathematicaでは，Tableの中身を{i,j}から{RandomReal[1], RandomReal[1]}に置き換えて最後の結果の/nを除く。

［１］地域内距離（腰塚武志）