Scikit-learn 青空を再現描画



#始めに
![sky.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/391380/657d0aad-1502-4863-306e-aec881dd9e83.png)

画像などで、良く電線などが写っていて、何とかしたいと思うことがあります。
上記画像で電線など余計なものを、機械学習などの応用で消すことはできないか？
と思い、この鉄塔の右側を雲一つない晴天に書き換えることをScikit-learnを使って
やってみたいと思います。
Scikit-learnを使って技術的なお遊びで、気楽にやっていますので、scoreを取ったりとか
正解率がどうのということは記事中でははぶきます。


#処理
##処理概要
処理としては下記のように行いました。

1. 画像読み込み
1. 学習用データを作成
1. Pipelineで回帰計算インスタンスをRGBそれぞれ定義し学習させる。
1. 回帰により推定値をR,G<Bそれぞれ取得
1. 推定値で青空を描き、マスク処理で画像合成


##回帰計算
空のR,G,B値を求めて、空を描き直します。
機械学習として、座標x, 座標yを入力し、R,G,Bを出力するように
したいですが求められる推定値は1つなので、R,G,B毎に1つづつ
用意します。

+ B ← f(座標x, 座標y)
+ G ← f(座標x, 座標y)
+ R ← f(座標x, 座標y)


**f()**は**Pipeline()**でまとめられた回帰計算処理です。

##サンプリング
サンプリングは下記位置23か所で行います。
![skymark.jpg](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/391380/936ff3d0-bd97-c57c-7f57-9349f2137639.jpeg)

**Pythonコード**

```python
mes_pnts = np.array([
    (174,  6),( 11, 12),(362, 24),(233, 28),( 66, 29),
    (207, 32),(  9, 54),( 96, 62),(342, 73),(289, 92),
    (375,112),(387,152),(172,155),(230,163),(191,168),
    (176,188),(345,195),(  9,222),(383,230),(291,236),
    (208,256),( 47,259),(156,262),])

# 機械学習用データを作成
rw,rh = 3, 3
x_r,x_g,x_b = [],[],[]
for ps in mes_pnts:
    x = int(ps[0]-1)
    y = int(ps[1]-1)

    # サンプリング座標を中心にして
    # 3x3の9画素の平均を取得する
    m = cv2.mean(img[y:y+rh, x:x+rw])
    x_r.append(m[2])    # R平均
    x_g.append(m[1])    # G平均
    x_b.append(m[0])    # B平均
```

サンプリング値は指定座標を中心に3×3の領域の平均値を使用しています。
画像処理で1画素で判断するのは危険です。


##Pipelineの使用
Pipelineは下記機能を組み合わせて使用します。

+ StandardScalerの標準化機能
+ PolynomialFeaturesの多項式機能
+ LinearRegressionの回帰分析機能

これを使うことで入力、学習の操作がいっぺんにできます。
上記１つ１つで１記事が書けるほどなので説明がしきれず、本記事の守備範囲を
超えてしまうと感じる為、詳細についてはインターネット、Qiita上に沢山あるので
詳しく知りたい方はググって見て下さい。

**参考**
&emsp;&emsp;[scikit-learn.org - sklearn.pipeline.Pipeline](https://scikit-learn.org/stable/modules/generated/sklearn.pipeline.Pipeline.html)
&emsp;&emsp;[scikit-learn.org - sklearn.preprocessing.PolynomialFeatures](https://scikit-learn.org/stable/modules/generated/sklearn.preprocessing.PolynomialFeatures.html?highlight=polynomialfeatures#sklearn.preprocessing.PolynomialFeatures)
&emsp;&emsp;[scikit-learn.org - sklearn.linear_model.LinearRegression](https://scikit-learn.org/stable/modules/generated/sklearn.linear_model.LinearRegression.html?highlight=linearregression#sklearn.linear_model.LinearRegression)
&emsp;&emsp;[scikit-learn.org - sklearn.preprocessing.StandardScaler](https://scikit-learn.org/stable/modules/generated/sklearn.preprocessing.StandardScaler.html?highlight=standardscaler#sklearn.preprocessing.StandardScaler)
<br>
**Rチャンネルを処理するPythonコード**

```python
from sklearn.pipeline import Pipeline
from sklearn.linear_model import LinearRegression 
from sklearn.preprocessing import PolynomialFeatures
from sklearn.preprocessing import StandardScaler

order = 5
pf_r = Pipeline([
        ("scaler", StandardScaler()),
        ("poly_features", PolynomialFeatures(degree=order,include_bias=False)),
        ("linear_reg", LinearRegression(fit_intercept=True)) ])
pf_r.fit(mes_pnts, x_r) # R学習

# 推定
prod_r = np.clip(pf_r.predict(wp),0,255).astype(np.int32)
```

##マスク処理
マスク画像としては下記画像を使います。白い部分(255)に青空を描きます。
![rchmask.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/391380/95f7719a-5458-9bd1-2bd6-70404811fded.png)

マスク処理の流れ

1. 反転マスク画像 ← NOT マスク画像
1. 合成用画像 ← 入力画像 AND 反転マスク画像
1. 結果画像 ← 生成青空画像 AND 合成用画像


#Pythonコード

```python
from sklearn.pipeline import Pipeline
from sklearn.linear_model import LinearRegression 
from sklearn.preprocessing import PolynomialFeatures
from sklearn.preprocessing import StandardScaler
import numpy as np
import math
import cv2

def main():
    # 鉄塔画像の読み込み
    img = cv2.imread('sky.png', cv2.IMREAD_COLOR)
    # マスク画像の読み込み
    mask = cv2.imread('rchmask.png', cv2.IMREAD_GRAYSCALE)
    mask = cv2.bitwise_not(mask)
    height, width = img.shape[:2]
    mes_pnts = np.array([   # サンプリング座標
        (174,  6),( 11, 12),(362, 24),(233, 28),( 66, 29),
        (207, 32),(  9, 54),( 96, 62),(342, 73),(289, 92),
        (375,112),(387,152),(172,155),(230,163),(191,168),
        (176,188),(345,195),(  9,222),(383,230),(291,236),
        (208,256),( 47,259),(156,262),])

    # 機械学習用データを作成
    rw,rh = 3, 3
    x_r,x_g,x_b = [],[],[]
    for ps in mes_pnts:
        x = int(ps[0]-1)
        y = int(ps[1]-1)
        # サンプリング座標を中心にして
        # 3x3の9画素の平均を取得する
        m = cv2.mean(img[y:y+rh, x:x+rw])
        x_r.append(m[2])    # R平均
        x_g.append(m[1])    # G平均
        x_b.append(m[0])    # B平均

    order = 5
    pf_r = Pipeline([
            ("scaler", StandardScaler()),
            ("poly_features", PolynomialFeatures(degree=order,include_bias=False)),
            ("linear_reg", LinearRegression(fit_intercept=True)) ])
    pf_r.fit(mes_pnts, x_r) # R学習
    pf_g = Pipeline([
            ("scaler", StandardScaler()),
            ("poly_features", PolynomialFeatures(degree=order,include_bias=False)),
            ("linear_reg", LinearRegression(fit_intercept=True)) ])
    pf_g.fit(mes_pnts, x_g) # G学習
    pf_b = Pipeline([
            ("scaler", StandardScaler()),
            ("poly_features", PolynomialFeatures(degree=order,include_bias=False)),
            ("linear_reg", LinearRegression(fit_intercept=True)) ])
    pf_b.fit(mes_pnts, x_b) # B学習

    # マスクより空を再描画する
    ya, xa = np.where(mask == 0)
    wp = np.hstack((xa.reshape(-1,1), ya.reshape(-1,1)))

    # 推定
    prod_r = np.clip(pf_r.predict(wp),0,255).astype(np.int32)
    prod_g = np.clip(pf_g.predict(wp),0,255).astype(np.int32)
    prod_b = np.clip(pf_b.predict(wp),0,255).astype(np.int32)

    # 結果画像
    result = cv2.bitwise_and(img, cv2.merge((mask,mask,mask)))
    result[ya, xa, 0] = prod_b
    result[ya, xa, 1] = prod_g
    result[ya, xa, 2] = prod_r

    # 表示
    cv2.imshow("original", img)
    cv2.imshow("result", result)
    cv2.waitKey(0)

if __name__ == '__main__':
    main()
```

sky.pngファイルが入力画像です。


#結果
上記コードを実行することで、下記画像が得られます。

![skyresult2.jpg](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/391380/4f2f7a42-671c-fa07-1f66-1aeff11b9498.jpeg)


#おわりに
画像としては、いいかんじにできたかなとは思います。際の部分ではイマイチ感ありますが。
この処理では、実は入力のサンプリング位置の確定や、マスク画像を作る方が
難しかったりします。その辺の方法についても、今後いろいろ試したいなとは思います。
お読み下さりありがとうございました。


#参考
[Emotion Explorer - 機械学習で空いっぱいの電線、鉄塔を消すことを試す](https://emotionexplorer.blog.fc2.com/blog-entry-330.html)
[Emotion Explorer - Scikit-learn Pipelineを試す](https://emotionexplorer.blog.fc2.com/blog-entry-284.html)
[Emotion Explorer - Scikit-learn 多項式回帰と交互作用項](https://emotionexplorer.blog.fc2.com/blog-entry-301.html)

#写真引用
鉄塔と青空 [photoAC](https://www.photo-ac.com)掲載 [ボッチ写真家](https://www.photo-ac.com/profile/2589747)さんの写真