質問をすることでしか得られない、回答やアドバイスがある。

15分調べてもわからないことは、質問しよう!

新規登録して質問してみよう
ただいま回答率
85.50%
Python 3.x

Python 3はPythonプログラミング言語の最新バージョンであり、2008年12月3日にリリースされました。

Q&A

解決済

8回答

3281閲覧

numpyの連続した数値列のカウント

p_pp

総合スコア17

Python 3.x

Python 3はPythonプログラミング言語の最新バージョンであり、2008年12月3日にリリースされました。

0グッド

4クリップ

投稿2018/04/12 13:04

編集2018/04/12 14:34

python初心者です
二値化した画像(.jpg)を配列に入れました
配列が

[1 1 1 1 1 1 1 1 0 1 0 0 0 0 0]
[1 1 1 1 1 1 1 1 0 0 0 0 0 0 0]
[1 1 1 1 1 1 1 1 1 1 0 0 0 0 0]

のようになっているとき,連続した1の数をカウントするソースコードを教えていただきたいです
1の数を数えるソースコードは調べれば出てくるのですが、連続した数字を数えるソースコードは出てきません

よろしくお願いいたします

実際には画像ファイルが多く,ピクセル数も多いため時間がかかってしまいます
1,0の境界の位置を探してそのピクセル数を得るのに時間短縮できる方法があれば,教えていただけると有難いです
最終的には写真から輪郭形状のグラフを作りたいと思っています

気になる質問をクリップする

クリップした質問は、後からいつでもMYページで確認できます。

またクリップした質問に回答があった際、通知やメールを受け取ることができます。

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

umyu

2018/04/12 13:19

最終的に行いたいことはなにですか?輪郭抽出とかでしょうか?
p_pp

2018/04/12 13:28

そうです!写真から最終的にはグラフの作成を考えています
umyu

2018/04/12 13:33

>itakovichさんへ 質問文にその旨を追加してくださいな。
p_pp

2018/04/12 14:33

分かりにくくて、すみません… 追加しておきます
guest

回答8

0

https://www.geeksforgeeks.org/maximum-consecutive-ones-or-zeros-in-a-binary-array/

普通に考えて全部のピクセルを見ないといけないのでO(N)になります。(worstケース)

最も簡単な方法は並列化です。
今回の場合、単純並列なので、立てることのできるスレッドの数分の1に計算時間を減らすことができます。


あとはアルゴリズムの話になりますが、例えばk番目まで見て、今のマスは0でこれまでの最長の1たちは100個でした。

そうすれば、とりあえず次は、k+100を見てみることにします。
こいつが1であればk+50を見てみることにします。
そして…

なんとなくお察しいただけるかと思いますが、このように探索することで枝切りの可能性が増えます。
例えばk+50の時点で0に出逢えば、k+1k+49k+51k+99はもう探索する必要はありません。

このようにすることで常に153620481000個のピクセルを見なければいけなかったところを、より少ないピクセルを探索すれば良いことになります。

後はデータの性質に依存します。
01010101010101010101010101010のように並んでいるのであれば判定が必要な分遅くなります。
11111111110000011100000000011のようであれば、大幅に探索を減らすことができます。


ブロック分けしてdivide-and-conquerでもよいですが。

↓ ↓↓  ↓ 0000001 0111011 ... 0000010 1111011 ... 1000010 0000000 ... 1100001 1000001 ... 1000010 0111110 ... ...

比較。(更新:groupby, numpy.where, backward, forward)

python

1from itertools import groupby 2import numpy as np 3import random 4 5from contextlib import contextmanager 6import time 7@contextmanager 8def timer(name): 9 t0 = time.time() 10 yield 11 print(f'[{name}] done in {time.time() - t0:.5f} s') 12 13slen = 1300 14tslen = 1500 15N = 10000 16 17ss = [] 18for _i in range(N): 19 random.seed(_i) 20 s_ = np.array([1]*(tslen-slen) + [random.randint(0, 1) for _ in range(slen)]) 21 ss.append(s_) 22 23with timer('groupby'): 24 mm0 = [] 25 for _i in range(N): 26 s = ss[_i] 27 ans = [len(list(g)) for k, g in groupby(s)] 28 m = max(ans[::2]) 29 mm0.append(m) 30 31with timer('numpy.where'): 32 mm1 = [] 33 for _i in range(N): 34 a = ss[_i] 35 offset = 1 36 diff_a = a[offset:] - a[:-offset] 37 38 pos1to0 = np.where(diff_a ==-1)[0] + offset 39 pos0to1 = np.where(diff_a ==+1)[0] + offset 40 41 pos0to1 = np.append([0], pos0to1) 42 m = (pos1to0 - pos0to1[:len(pos1to0)]).max() 43 mm1.append(m) 44 45with timer('backward-search'): 46 def forward_search(s): 47 j = 0 48 while j < len(s) and s[j] == 1: 49 j += 1 50 return j 51 52 def backward_search(s, m): 53 if m >= len(s): 54 return False, len(s), 0 55 t = s[m::-1] 56 c = forward_search(t) 57 return True, m, c 58 59 mm2 = [] 60 for _i in range(N): 61 s = ss[_i] 62 m = 0 63 i = forward_search(s) 64 m = i 65 while i < len(s): 66 fforward, j, c = backward_search(s[i:], m) 67 i += j 68 if fforward: 69 j = forward_search(s[i:]) 70 c += j 71 i += j 72 if c > m: 73 m = c 74 mm2.append(m) 75 76with timer('forward-search'): 77 mm3 = [] 78 for _i in range(N): 79 s = ss[_i] 80 m = 0 81 c = 0 82 for v in s: 83 if v == 1: 84 c += 1 85 else: 86 if c > m: 87 m = c 88 c = 0 89 if c > m: 90 m = c 91 mm3.append(m) 92 93assert np.all(mm0 == mm1) 94assert np.all(mm0 == mm2) 95assert np.all(mm0 == mm3)

雑談

高速化を行う際には、

  1. アルゴリズムの改良
  2. コーディングの改良

に分けて考えるとスッキリします。

今回の場合、最低でもO(N)で各ピクセルを見ないと答えが出ないような問題です。
すると、探索の打ち切りの可能性を検討した後、コーディングの改良を検討します。

(目的とあっていない可能性がありますが、
最長の1の塊を探索するのであれば、最長になりえない領域をできるだけスキップする作戦を考えます。)

次に書き方を工夫します。
(ブロック化してSIMD命令を利用するのが、C++などのコンパイル言語を使う際に検討する必要があります。)
Pythonの場合、numpy、opencvのライブラリを出来る限り利用できるようにするのが近道です。
すると、numpyのメソッド、opencvのメソッド、pythonのbuilt-inメソッドを使用することを検討すべきです。

一個ずつピクセルを見ていくより、ベクトルとみなしてsumが速いことは予想できます。
ならば、sumだけで済むように前処理を検討するのが良さそうです。

opencvのオープニングを使うことで先にノイズを消してから、sumを行うのが一番速そうです。

他にもバイナリとして扱ってビット演算を用いることで高速化することも考えられます。
Pythonでそれをやるのは少しちぐはぐ感が残ってしまいます。
C++などで書き換えるだけで、もとのアルゴリズムで高速に処理できることが見込まれるためです。

(重複するようなケースが出現する場合、キャッシュをうまく使うことも重要です。
https://pypi.python.org/pypi/fastcache/0.4.3
メモリと演算のバランスをチューニングするがよいコード書くために大切です。)

高速化することも大事ですが、楽して高速化することも大事です。

投稿2018/04/12 13:40

編集2018/04/15 13:41
mkgrei

総合スコア8560

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

p_pp

2018/04/12 14:43

分かりやすい説明ありがとうございます! さっそく実装してみます
p_pp

2018/04/12 19:21

やってみました! 意外と難しくて、なかなかうまくいきません データの性質では 111111111111111111111110000000000000000000 111111111111111111111111100000100000000000 111111111111111000000000001000000000000000 のように一列目はすべて1から始まっています ノイズをカウントしたくないので連続しているという1の数を求めたいという標記にしました よろしければ探索を減らすソースコードを乗せていただけないでしょうか
mkgrei

2018/04/15 05:27

今気づいたのですが、 連続したブロックのうち最長のものがほしいのか、 1個以上連続で出現している1の総数がほしいのか、 どちらでしょう?
mkgrei

2018/04/15 05:43

前者ならデータの性質の応じて最後に追記したコードで処理できます。 連続する1が少ない場合、素直に順番に探索するのが最速です。 連続する1が多い場合、枝切りをすることで速めることができます。 ノイズをカウントせずに、1個以上連続で出現している1の総数がほしいのであれば、 http://labs.eecs.tottori-u.ac.jp/sd/Member/oyamada/OpenCV/html/py_tutorials/py_imgproc/py_morphological_ops/py_morphological_ops.html オープニングを行うことでノイズを消せます。 後に、和を取ることによって1の数を得ることができます。
p_pp

2018/04/16 07:11

ありがとうございます!! こっちもやってみます!!
guest

0

ベストアンサー

0と1の境界線が欲しいなら、変化点を探して距離を測ればいいんじゃないかなぁ。変化点なんで、まず微分します。

python

1a = np.array([1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 0, 0, 0, 0]) 2offset = 1 3diff_a = a[offset:] - a[:-offset]

1->0は-1の変化、0->1は+1の変化なので、それぞれの場所をnp.argwhereを使って抽出します。

pos1to0 = np.where(diff_a ==-1)[0] + offset # -> [8, 10] pos0to1 = np.where(diff_a ==+1)[0] + offset # -> [9]

右端が必ず0で終わると仮定すれば、pos0to1を微調整してインデックスの差分を取れば1が連続する長さを測れます。

pos0to1 = np.append([0], pos0to1) len11 = pos1to0 - pos0to1

pythonのfor-loopが遅いので行列のインデックスを駆使した方法です。理解してコーディングしないと収拾がつかなくなるので、適度な用法用量を守ってご参考にしてください。


追記

mkgreyさんのコードそのまま回したケース
[groupby] done in 0.39933 s
[forward-search] done in 0.12787 s
[backward-search] done in 0.05486 s

最大の長さを返すよう若干修正したバージョン

In [51]: slen = 1000 ...: tslen = 1200 ...: N = 1000 ...: ...: random.seed(0) ...: s_org = [1]*(tslen-slen) + [random.randint(0, 1) for _ in range(1000)] ...: In [52]: with timer("numpy"): ...: a = np.array(s_org) ...: offset = 1 ...: diff_a = a[1:] - a[:-1] ...: pos1to0 = np.where(diff_a ==-1)[0] + offset ...: pos0to1 = np.where(diff_a ==+1)[0] + offset ...: pos0to1 = np.append([0], pos0to1) ...: len11 = (pos1to0 - pos0to1).max() ...: [numpy] done in 0.00026 s In [53]: assert len11 == 202

投稿2018/04/15 08:37

編集2018/04/15 12:17
tachikoma

総合スコア3601

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

mkgrei

2018/04/15 09:21

文句をつけるようで申し訳ないのですが、手元の環境でやってみたところfor文を回したほうが速かったです。 私の実装に自信がないので、ご確認していただいてよろしいでしょうか。
tachikoma

2018/04/15 09:58

なるへそ。コメントありがとうございます。しばしお待ちを。
tachikoma

2018/04/15 12:18

mkgreiさん、比較してみました。同じものを比べてますかね??np.arrayを前もって変換しておけば実際の計算時間はこれの半分程度で済みました。
mkgrei

2018/04/15 12:32

ご確認していただいてありがとうございます。 a = np.array([1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 0, 0, 0, 0]) で比較していました。 numpyの方が長さが増えても計算時間がほとんど変わらず、速いですね。
tachikoma

2018/04/15 12:57

numpyのテクニックですね(ちなみにMATLABも同様)。アルゴリズムを工夫しても実装であしを引っ張られてしまうのは辛いところ。
mkgrei

2018/04/15 13:15

whereを使えば、最長だけでなくすべての区間も取り出せるのがいいですね。 lookupなら素直に組んでもワンチャンあるかとも思いましたが、ダメですね。
p_pp

2018/04/16 06:09

返信遅れてすいません 分かりやすい解説ありがとうございます!!
p_pp

2018/04/19 09:23

np.whereで-1を探しても出てこないのですが、なぜでしょうか
tachikoma

2018/04/19 09:35

値が変化してない・・・ということは無いですか?
p_pp

2018/04/20 02:58

確認してみましたが値は変化していました…
p_pp

2018/04/20 07:01

from PIL import Image import numpy as np import matplotlib.pyplot as plt np.set_printoptions(threshold=np.inf) #閾値 t=100 for i in range(0,1): filename='C://toisi//f_%05d.jpg'%(i) im = np.array(Image.open(filename).convert('L')) #二値化 im[im<t]=0 im[im>=t]=1 #print(im.shape, im.dtype) #print(im[1,1250]) #0の数をカウント for i in range(1535): a = im[i] offset = 1 diff_a = a[offset:] - a[:-offset] pos0to1 = np.where(diff_a ==+1)[0] + offset pos1to0 = np.where(diff_a ==-1)[0] + offset #print(pos0to1) #print(pos1to0) print(diff_a)
p_pp

2018/04/20 07:06

diff_aをプリントしてみてみると [0,0,0,0,0,0,0,0,0,0,0,0,255,1,1,0,0,0,0,0,0,0] のような感じになってしまい、出るはずのない255という値がでてきてしまいました
LouiS0616

2018/04/20 07:22

dtypeがnp.uint8なのでしょう。アンダーフローが発生しています。 この場合気にせず diff_a == 255 の箇所を利用すれば良いかと思います。
tachikoma

2018/04/20 07:31

フォローありがとうございます。画像を読み込んだ結果がuint8になっているためで、その型uint8同士の演算もuint8になって、結果アンダーフローして255に飛んだわけですか。im=im.astype(np.int)でint型にするか、LouiS0616さんのコメントにあるとおり255を利用するかですね。迂闊でした。
p_pp

2018/04/20 09:03

ありがとうございます!! できました! 少し聞きたいのですがpos0to1のときは[9]でpos1to0の中身が[]のように何もないときpos1to0-pos0to1でエラーが表示されてしまうのですがif文で無理やり入れてしまうなどの対処法でよかったでしょうか 度々すいません。よろしくお願いいたします
mkgrei

2018/04/20 09:40

m = (pos1to0 - pos0to1[:len(pos1to0)]).max() どっちが短いか見極めて、スライスしてください。
p_pp

2018/04/25 11:46

ありがとうございます!!
guest

0

一列目はすべて1から始まっています

本来の目的から考えると、私なら左端を塗りつぶした画像を新たに作成します。
OpenCVを使えば内部はネイティブ処理してくれると思うので速度も期待できます。
以下例では左端の黒を白ですが、逆も同様の考えでできます。

Python

1import cv2 2img = cv2.imread( 'bin.bmp', cv2.IMREAD_GRAYSCALE) # 簡単のため2値ビットマップを読込 3print(img.shape) 4# y方向に左端が黒を白に塗りつぶし 5for r in range(img.shape[0]): 6 if img[r,0] == 0: 7 print('fill',r) 8 cv2.floodFill(img,None,(0,r),255) 9cv2.imwrite('out.jpg',img)

bin.bmp
イメージ説明
out.jpg
イメージ説明

投稿2018/04/13 03:48

can110

総合スコア38233

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

0

初めて1が来たところから0が来るまでfor文やwhile文でカウントし、カウントした数を別のリストに格納する。if文で評価し、0だったらbreak文で切る

これでどうでしょうか?

投稿2018/04/12 13:13

退会済みユーザー

退会済みユーザー

総合スコア0

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

p_pp

2018/04/12 13:17

早速の回答ありがとうございます! 一度それでやったのですが時間がかかってしまいました 1536×2048の画像で,1000個のファイル処理をする予定です
退会済みユーザー

退会済みユーザー

2018/04/12 13:24

なるほど... もっと効率的な方法を考えておきます。
p_pp

2018/04/12 13:29

すみません ありがとうございます!!
退会済みユーザー

退会済みユーザー

2018/04/12 13:48 編集

今思いついたアルゴリズムを書きたいと思います。 同じ数字が続くということはその二つの整数の差は"0"です。0~1または1~0が変わるということは 差の絶対値は1です。この差の絶対値が1というところが変わる境目なのでこれを数えるのはどうでしょうか? 差が1のところと-1のところを別々にカウントするという意味です。 (データが多いのでこれも時間がかかるかもしれません)
p_pp

2018/04/12 14:37

ありがとうございます!! やってみます!
p_pp

2018/04/12 19:16

やってみましたがあまり変わらなかったです… ありがとうございます!
退会済みユーザー

退会済みユーザー

2018/04/13 00:52

すみません、もう少し考えたらカウントではなく差が-1や1になるところの"リストの番号"をとりだす だけでいいと思いました。
退会済みユーザー

退会済みユーザー

2018/04/13 01:04

あと、具体的な方法はまだですが別方法で"0と1"のみなので「2進数を使う」も検討中です。
guest

0

RLE (Run Length Encoding)を高速にできれば、境界をグラフ化するのも簡単になりそうですね。

Pure PythonでRLEを書くにはitertools.groupbyを使うのが正攻法らしいです。Pure Pythonに限らなければもっと速く書くことはできると思います。

python

1from itertools import groupby 2[(k, len(list(g))) for k, g in groupby("111100010101011111000")] 3""" 4[('1', 4), 5 ('0', 3), 6 ('1', 1), 7 ('0', 1), 8 ('1', 1), 9 ('0', 1), 10 ('1', 1), 11 ('0', 1), 12 ('1', 5), 13 ('0', 3)] 14"""

長さ1(非連続)の1を取り除くのはRLE後でも可能ですが、手間がかかるのでできるだけRLE前に別途画像処理系のアルゴリズムで取り除いておくほうが良いと思います。良くは知りませんがノイズ除去処理あたりだと思います。

投稿2018/04/15 02:38

YouheiSakurai

総合スコア6142

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

0

もしかしてnumpy配列を直接pythonから処理していますか?
それなら、cythonの力で高速化しちゃえば同じアルゴリズムで100倍くらいは期待できるはずです。

参考ページ:
NumPyとCythonを組み合わせると爆速! - Kesinの知見置き場

投稿2018/04/13 01:03

hayataka2049

総合スコア30933

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

0

シンプルに書けば速くなる‥なんて甘くはないだろうけど。

python

1import numpy as np 2 3pixels = np.array([[1,1,1,0,0,0],[1,0,0,0,0,0],[1,1,1,1,1,0],[1,1,1,1,1,1]]) 4 5def search0(line): 6 try: 7 return line.index(0) 8 except: 9 return 0 10 11print([search0(list(line)) for line in pixels]) 12# [3, 1, 5, 0]

絶対に0が含まれているなら、

python

1pixels = np.array([[1,1,1,0,0,0],[1,0,0,0,0,0],[1,1,1,1,1,0]]) 2print([list(line).index(0) for line in pixels]) 3# [3, 1, 50]

投稿2018/04/13 01:47

fuzzball

総合スコア16731

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

0

元の配列からひとつずらした配列を作成して、ベクトル同士の引き算をします。結果のベクトルの要素がゼロならひとつ前の値と連続状態にあることになります。
処理はnumpyのarrayにすれば簡単かと思います。2次元ベクトルと見なしてしまえばこの差の計算は1回で終わります。
arrayのsliceを使ったほうが簡単かもしれません

投稿2018/04/12 21:11

編集2018/04/12 21:25
R.Shigemori

総合スコア3376

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

あなたの回答

tips

太字

斜体

打ち消し線

見出し

引用テキストの挿入

コードの挿入

リンクの挿入

リストの挿入

番号リストの挿入

表の挿入

水平線の挿入

プレビュー

15分調べてもわからないことは
teratailで質問しよう!

ただいまの回答率
85.50%

質問をまとめることで
思考を整理して素早く解決

テンプレート機能で
簡単に質問をまとめる

質問する

関連した質問