2016-06-23 10 views
5

私は、numpyを使用して最適化したい、好ましくはループを削除したい次のコードを用意しました。私はそれにアプローチする方法を見ることができないので、どんな示唆も役に立つでしょう。ループの最適化/削除

インデックスは、(N、2)個数の整数配列であり、Nは数百万にすることができます。コードが行うことは、最初の列に繰り返しインデックスを見つけることです。これらのインデックスについては、2番目の列に対応する2つのインデックスのすべての組み合わせを作成します。その後、最初の列のインデックスと一緒にそれらを収集します。

index_sets = [] 
uniques, counts = np.unique(indices[:,0], return_counts=True) 
potentials = uniques[counts > 1] 
for p in potentials: 
    correspondents = indices[(indices[:,0] == p),1] 
    combs = np.vstack(list(combinations(correspondents, 2))) 
    combs = np.hstack((np.tile(p, (combs.shape[0], 1)), combs)) 
    index_sets.append(combs) 
+0

ネットワークに問題がありますので、おそらく 'networkx'モジュールを見てください。 – Divakar

答えて

1

ここでは、Nに対してベクトル化された解法があります。まだforループが1つ含まれていることに注意してください。ただし、それぞれの 'キーの多重度グループ'にはループがあります。多くても数ダース)。

N = 1.000.000の場合、ランタイムはPC上で1秒の大きさです。

import numpy_indexed as npi 
N = 1000000 
indices = np.random.randint(0, N/10, size=(N, 2)) 

def combinations(x): 
    """vectorized computation of combinations for an array of sequences of equal length 

    Parameters 
    ---------- 
    x : ndarray, [..., n_items] 

    Returns 
    ------- 
    ndarray, [..., n_items * (n_items - 1)/2, 2] 
    """ 
    return np.rollaxis(x[..., np.triu_indices(x.shape[-1], 1)], -2, x.ndim+1) 

def process(indices): 
    """process a subgroup of indices, all having equal multiplicity 

    Parameters 
    ---------- 
    indices : ndarray, [n, 2] 

    Returns 
    ------- 
    ndarray, [m, 3] 
    """ 
    keys, vals = npi.group_by(indices[:, 0], indices[:, 1]) 
    combs = combinations(vals) 
    keys = np.repeat(keys, combs.shape[1]) 
    return np.concatenate([keys[:, None], combs.reshape(-1, 2)], axis=1) 

index_groups = npi.group_by(npi.multiplicity(indices[:, 0])).split(indices) 
result = np.concatenate([process(ind) for ind in index_groups]) 

免責事項:私はnumpy_indexedパッケージの作者です。

+0

ありがとう、私はこの解決策の性能をテストする必要があります。 – martinako

+0

試してみましたか?それが実際にどのようにうまくいったか聞いてみるのは興味深い。 –

+0

申し訳ありませんまだ試してみることができませんでした。より優先度の高い別のタスクに切り替える必要がありましたが、このコードを使用しています。私はその日の終わりにそれをテストし、報告しようとします。 – martinako

2

ほとんど改善が示唆され得る:

  • は、我々は、各グループに対応する組合せを記憶するために必要な行の推定数を事前に計算することができるため、出力アレイを初期化します。 N要素では、可能な組み合わせの合計数はN*(N-1)/2となり、各グループの組み合わせの長さがわかります。さらに、出力配列の行の総数は、これらの間隔の長さの合計になります。

  • ループに入る前にできるだけ多くのものをベクトル化して事前計算します。

  • ループを使用して、不規則なパターンをベクトル化できないため、組み合わせを取得します。タイリングをシミュレートするにはnp.repeatを使用し、beforeループを実行して各グループの最初の要素、つまり出力配列の最初の列を取得します。

だから、心の中ですべてのものを改良して、実装は次のようになります -

# Remove rows with counts == 1 
_,idx, counts = np.unique(indices[:,0], return_index=True, return_counts=True) 
indices = np.delete(indices,idx[counts==1],axis=0) 

# Decide the starting indices of corresponding to start of new groups 
# charaterized by new elements along the sorted first column 
start_idx = np.unique(indices[:,0], return_index=True)[1] 
all_idx = np.append(start_idx,indices.shape[0]) 

# Get interval lengths that are required to store pairwise combinations 
# of each group for unique ID from column-0 
interval_lens = np.array([item*(item-1)/2 for item in np.diff(all_idx)]) 

# Setup output array and set the first column as a repeated array 
out = np.zeros((interval_lens.sum(),3),dtype=int) 
out[:,0] = np.repeat(indices[start_idx,0],interval_lens) 

# Decide the start-stop indices for storing into output array 
ssidx = np.append(0,np.cumsum(interval_lens)) 

# Finally run a loop gto store all the combinations into initialized o/p array 
for i in range(idx.size): 
    out[ssidx[i]:ssidx[i+1],1:] = \ 
    np.vstack(combinations(indices[all_idx[i]:all_idx[i+1],1],2)) 

を出力配列は大きな(M, 3)形状の配列でも、配列のリストに分割されませんのでご注意ください元のコードによって生成されます。それでも必要な場合は、同じものとしてnp.splitを使用できます。

また、クイックランタイムテストでは、提案されたコードではあまり改善されていないことが示唆されています。したがって、おそらくランタイムの大部分は組み合わせを得るために費やされます。したがって、そのような接続ベースの問題に特に適しているnetworkxの代わりの方法が適しているようです。

+0

私はこれを私の答えと比較してベンチマークを試みましたが、N = 1.000.000のMemoryErrorを返しました:) –