Pythonで配列内の要素をを分割したい

Question

前提・実現したいこと

以下の様なアルゴリズムを実装したいです。

1.配列に01までのランダムな小数が入っている。
2.配列から要素を1つ持ってくる。
3.その要素が01をN等分したどこに入るかを分配する。
(例:取り出した要素=0.598,等分数=100の時、0.59<=要素<=0.6なので59番目)
4.新たな配列を用意し、3で決まった番号の回数を1個増やす
(例:先程0.598は59番目に入るのが分かったので、新たな配列の59番目を＋1する。)
5.24を最初の要素数行う。
6.4,5で作成した配列の要素を持ってくる。
7.その要素数分だけ等分し、どこに入るかを考える。
(例:59番目の配列には5個の要素があるとすると、0.590.6の間を5等分し、0.58<=0.598<=0.6より、0.598はその中の5番目に入る)
8.7を分配後の要素の数だけ繰り返し、再等分された中に入っているかどうか割合を考える。それを更に新しい配列に順次追加する。
(例:5等分され、1,3,5番目に数が1個以上あり、2,4番目に入っていない場合,3/5となる。)
9.6~8を元々のN等分行う。

6以降の追記です。
5までで、59番目に入った配列の要素が5個、その内容を[0.5910,0.5930,0.5930,0.5960,0.5980]とします。
ここから、0.590.60に入っているこの5個の要素を更に細かく見ていくこととなります。
再等分数は、要素が5個の為0.590.60を5等分します。この5等分に先程の要素が入っているかを考えます。
0.5910=1番目
0.5930=2番目
0.5930=2番目
0.5960=3番目
0.5980=4番目
すると、各区間に入っている数は

区間	要素が入っているかどうか
0.5900~0.5920	○
0.5921~0.5940	○
0.5941~0.5960	○
0.5961~0.5980	○
0.5981~0.6000	×

となります。ここから、5等分された区間の中にどの区間数が入ってるかを考えます。
今回の場合、5区間中に4区間は要素が存在しているので、最終的な出力は4/5となります。
この4/5を新しい配列の59番目に挿入します。

試したこと

5までは自力でやってみましたが、6からが分かりません。
1の配列作成は省略しています。(orbitに入っています。)

import numpy as np
import math

   N=100 #N:等分数。

   for i in range(len(orbit)):
        temporary = 0
        temporary = orbit[i]
        temporary = abs(math.floor(temporary*N))
        pdcount[temporary] = pdcount[temporary]+1

3の補足ですが、取り出した点をN倍し、小数を切り落としすればそのまま分配先の番号になることを利用しています。
(例:0.598,N=100の時、0.598*100=59.8,つまり59となり、pdcount[59]の要素を+1している。)

補足情報（FW/ツールのバージョンなど）

Python 3.94

Comments

[melian]

Numpy はお使いにならないのでしょうか？

[ntsk2543]

Nmupyは使用しています。追記しておきます。

[退会済みユーザー]

最終的にどのような結果が欲しいのか分かりません。ダミーサンプルでかまわないので0.59～0.6の間に10個くらい値があるケースで説明してもらえませんか?

[ntsk2543]

分かりました。追記しておきます。

[melian]

5 までの操作はいわゆるヒストグラムかと思いますので、

pdcount = np.histogram(orbit, bins=N, range=(0.0, 1.0))[0]

かと。

[退会済みユーザー]

更新ありがとうございます。重要な部分なので細かい点を確認させてください。7の区間を分割する部分は 0.590以上～0.592未満、0.592以上～0.594未満… みたいな感じじゃなくていいんですか?現在の例だと 0.59201などがどこにも含まれないことになります。

[ntsk2543]

区間の分割に関してはftlobwさんの言う通りです、あいまいな表現、失礼いたしました。

melianさんの質問ですが、記述の通りでした。ありがとうございました。

Answer 1

こんなことでしょうか？

import numpy as np
import random

random.seed(110)
orbit = np.asarray([random.uniform(0,1) for _ in range(10000000)])
print('start')

N = 100 # 分割数
hist, bins = np.histogram(orbit, range=(0, 1), bins=N) 
print(hist)
print(bins)

results = []

for i in range(len(bins)-1):
  r_min = bins[i]
  r_max = bins[i+1]

  pick_elems = orbit[(r_min <= orbit) & (orbit < r_max)]

  # print(r_min, r_max)
  # print(pick_elems)

  m_hist, _ = np.histogram(pick_elems, range=(r_min, r_max), bins=len(pick_elems)) 
  # print(m_hist)

  m_p = np.count_nonzero(m_hist)/len(m_hist)

  results.append(m_p)

print(len(results))
print(results)

[補足]
少しpythonっぽくしてみました

 !pip install more-itertools

import numpy as np
import random
import more_itertools

random.seed(110)
orbit = np.asarray([random.uniform(0,1) for _ in range(10000000)])
print('start')

N = 100

hist, bins = np.histogram(orbit, range=(0, 1), bins=N) 

def get_ratio(r_min:float, r_max:float, arr:np.ndarray)->float:
  pick_elems = arr[(r_min <= arr) & (arr < r_max)]
  m_hist, _ = np.histogram(pick_elems, range=(r_min, r_max), bins=len(pick_elems)) 
  return np.count_nonzero(m_hist)/len(m_hist)

results = [get_ratio(r[0], r[1], orbit) for r in more_itertools.windowed(bins,2,step = 1)]

print(len(results))
print(results)

Answer 2

import numpy as np
from numpy.random import default_rng

M = int(1e+7)
N = 100
num_range = (0.0, 1.0)

rg = default_rng()
orbit = rg.uniform(*num_range, M)

orbit_org = orbit.copy()
orbit.sort()
pdcount, bin_edge = np.histogram(orbit, bins=N, range=num_range)
bins = np.split(orbit, pdcount.cumsum()[:-1])

ratio = np.zeros(N)
for i, data in enumerate(bins):
  hist, _ = np.histogram(data, bins=pdcount[i], range=(bin_edge[i], bin_edge[i+1]))
  ratio[i] = hist[hist>0].size/pdcount[i]

print(ratio)
#
[0.63314754 0.63316065 0.6325016  0.63304742 0.6318326  0.63409262
 0.6309113  0.6328047  0.63185449 0.63081087 0.63180382 0.63368009
 0.63418594 0.63202644 0.630325   0.62935527 0.63230635 0.63173384
 0.63192124 0.63125318 0.63181714 0.63249917 0.63271949 0.63248817
 0.63252976 0.63157    0.63182343 0.63087606 0.63256839 0.6310312
                                 :
                                 :

以下の環境での実行時間は N = 1e+7 で 0.66±0.02 秒程度になりました(timeit モジュールを利用して計測)。

$ lscpu | grep -E '^(Architecture|Model name)'
Architecture: x86_64
  Model name: Intel(R) Core(TM) i5-8500T CPU @ 2.10GHz

$ lsb_release -ir
Distributor ID: Ubuntu
Release:        21.04

$ python3 --version
Python 3.9.5

$ python3 -c 'import numpy;print(numpy.__version__)'
1.21.4

Answer 3

処理を共通化してみました。

import numpy as np

N = 15

def classify(amin, amax, a, n):
    linspace = np.linspace(amin, amax, n+1)
    return [(pmin, pmax, orbit[(orbit >= pmin) & (orbit <pmax)]) for pmin, pmax in zip(linspace, linspace[1:])]

def ratio(amin, amax, a):
    alist = np.array([len(a) for pmin, pmax, a in classify(amin, amax, a, len(a))])
    return len(alist[alist>0])/len(alist)

result = [ratio(pmin, pmax, a) for pmin, pmax, a in classify(0, 1, orbit, N)]

実行例

>>> orbit = np.random.random(1000)
>>> result = [ratio(pmin, pmax, a) for pmin, pmax, a in classify(0, 1, orbit, N)]
>>> print(result)
[0.704225352112676, 0.6461538461538462, 0.6610169491525424, 0.6911764705882353, 0.6764705882352942, 0.6590909090909091, 0.6712328767123288, 0.5970149253731343, 0.5571428571428572, 0.5753424657534246, 0.6176470588235294, 0.6142857142857143, 0.59375, 0.6125, 0.6333333333333333]

Answer 4

記載された手順をそのまま素直に組んでみました。答え合わせはしていません。
要素数1000万で一様分布、N=100～10000でも10秒程度で終わります。
C++などでやるともっとはやくなると思います。

import math
import random

random.seed(110)
orbit = [random.uniform(0,1) for _ in range(10000000)]
print('start')
#print(orbit)

N = 100

# 各要素をクラス分け
slots = {}
for i, v in enumerate(orbit):
    n = math.floor(v*N)
    if n >= N:
        n = N-1
    if n not in slots:
        slots[n] = []
    slots[n].append(i) # 要素位置を保持

# 各クラスの計算
results = [0.0 for _ in range(N)] # 要素のないクラスの初期値は０
for n, lst in slots.items():
    exists = set()
    SUB_N = len(lst) # クラス内の要素数
    st = n / N # この区間の左端値
    for i in lst:
        v = (orbit[i]-st) * N # クラス内の値を0...1に正規化
        sub_n = math.floor(v*SUB_N)
        if sub_n >= SUB_N:
            sub_n = SUB_N-1
        exists.add(sub_n)

    results[n] = len(exists)/SUB_N

print(results)

Comments

[can110]

あ。クラス内の計算がおかしかった。これで合ってるかな？

Answer 5

0～1 を N等分するんだから 配列: slot[N] を用意し、
与えられた各要素:x に対し i/N <= x < (i+1)/N なら x を slot[i] に追加する。

...でいいんですよね?

Comments

[ntsk2543]

それが5までの内容になります。(質問がわかりずらくてすみません。)
ただ、本来ならこの後2次元配列のバージョンもあり、pythonで回すと処理に時間がかかってしまいます。(明記するのを忘れていましたが元々の配列の要素数は10万～1000万程ある。)
その為、各要素xを1回で振り分けることが可能になるようなアルゴリズムにしています。