はじめに

未来電子テクノロジーでインターンをしています、Sotaです。
プログラミング初心者のため、内容に誤りがあるかもしれません。
もし、誤りがあれば修正するので、どんどん指摘してください。

自分はc言語からプログラミングを勉強したので、%dや%fの表現に慣れています。
そのため、桁数指定などで、ストレスを感じたので、pythonでも%dや%fを使える表現をメモします。

pythonで%dや%sを使う

pythonでprintf文のような表現を使うには、
'文字列' % 変数 とすれば記述できます。
(例)

y = 10000
d = 100

print('100$は%d円' % y) #100$は10000万円
print('%d円は%d$' % (y, d)) #10000万円は100$

このような感じです。
もちろん文字列やfloat、桁数指定%.2fなども可能です。

format文

文字列メソッドformat()をでも似たような書き方ができます。
文字列中に置換フィールド{}が指定した引数と対応して置換されます。
文字列.format(変数)と書きます。
(例)

y = 10000
d = 100

print('100$は{}万円'.format(y)) #100$は10000万円
print('{}万円は{}$'.format(y, d)) #10000万円は100$

はじめに

未来電子テクノロジーでインターンをしています、Sotaです。
プログラミング初心者のため、内容に誤りがあるかもしれません。
もし、誤りがあれば修正するので、どんどん指摘してください。

自分はC言語からプログラミングを勉強したので、%dや%fの表現に慣れています。
そのため、桁数指定などで、ストレスを感じたので、Pythonでも%dや%fを使える表現をメモします。

Pythonで%dや%sを使う

Pythonでprintf関数のような表現を使うには、
'文字列' % 変数 とすれば記述できます。
(例)

y = 10000
d = 100

print('100$は%d円' % y) #100$は10000万円
print('%d円は%d$' % (y, d)) #10000万円は100$

このような感じです。
もちろん文字列やfloat、桁数指定%.2fなども可能です。

formatメソッド

文字列メソッドformat()でも似たような書き方ができます。
文字列中に置換フィールド{}が指定した引数と対応して置換されます。
文字列.format(変数)と書きます。
(例)

y = 10000
d = 100

print('100$は{}万円'.format(y)) #100$は10000万円
print('{}万円は{}$'.format(y, d)) #10000万円は100$

（ブログ記事からの転載）

[https://leetcode.com/problems/length-of-last-word/]

Given a string s consists of upper/lower-case alphabets and empty space characters ' ', return the length of last word in the string.

If the last word does not exist, return 0.

Note: A word is defined as a character sequence consists of non-space characters only.

Example:

Input: "Hello World"
Output: 5

Discussionで'This problem is not fun at all'とか言われちゃってましたが、、、LeetCodeは個人的にはなかなか難易度高い問題も多いので、こういう箸休め的問題は欲しいなと思います。
英単語で構成された文の中で、末尾の単語の長さを返す問題です。

解答・解説

解法1

リストlsに対して、2つのポインタとなる変数slow, fastを用意し、slowは末尾の単語の末尾を、fastは末尾の単語の先頭のインデックスを返すように計算します。
末尾の単語の長さを判定するので、リストを逆向きにループさせていきます。
まず、slowは初期値を-1とします。list[-1]とするとlistの末尾の要素を示すのは自明ですね。そこから逆向きにループさせて、ブランクでない位置まで辿ります。
通常、'Hellow World'のように末尾の単語の後ろにブランクはなく、slowは初期値の-1のまま（ループは進まない）になることが多いと思います。

ls = len(ls)
slow = -1
while slow >= -ls and s[slow] == ' ':
    slow -= 1

次に、fastは初期値をslowと同値とします。そこから逆向きにループさせて、ブランクの位置まで辿ります。

fast = slow
while fast >= -ls and s[fast] != ' ':
    fast -= 1

最後に、slow - fastとすれば、末尾の単語の長さを計算できます。
まとめると以下のようになります。

class Solution:
    def lengthOfLastWord(self, s: str) -> int:
        ls = len(s)
        slow = -1
        while slow >= -ls and s[slow] == ' ':
            slow -= 1
        fast = slow
        while fast >= -ls and s[fast] != ' ':
            fast -= 1
        return slow - fast

解法2

1 line solutionを示して、今日は終わります。rstripで文末の余計なブランクを剥ぎ取ってしまうという力技。

class Solution:
    def lengthOfLastWord(self, s):
        return len(s.rstrip(' ').split(' ')[-1])

動機

研究室内の課題で水汲み問題を解く問題が出たので始めにPyhtonで書き、Javascriptの勉強中だったのでJavascriptでも書いてWebサイト風にアレンジしてみました。これをサーバーで公開するのもなぁと思ったのでコピペしてブラウザで開いてみてください。

Pythonのコード

warer_jug.py

import queue

class Bottle:
    def __init__(self,one,two):
        self.one = 0 #一個目のボトル
        self.two = 0 #二個目のボトル
        self.one_max = one
        self.two_max = two

    def put(self,now_one,now_two,num,past_list): #oneのボトルを水で満たす
        copy = past_list.copy()
        if num ==1:
            now_one = self.one_max 
        elif num ==2:
            now_two = self.two_max
        for i in range(len(copy)):
            if copy[i] != [now_one,now_two]:
                if (i== len(copy)-1):
                    copy.append([now_one,now_two])
                    queue.put(copy)
                    break
            else:
                break

    def pour(self,now_one,now_two,num,past_list): #numから一方のボトルに水を移す
        copy = past_list.copy()
        if num==1: #oneからtwo
            now_two += now_one
            if now_two >= self.two_max: #あふれる場合,全部は入れない場合
                now_one = now_two - self.two_max
                now_two = self.two_max
            else:
                now_one = 0
        elif num ==2: #twoからone
            now_one += now_two
            if now_one >= self.one_max: #あふれる場合,全部は入れない場合
                now_two = now_one - self.one_max
                now_one = self.one_max
            else:
                now_two = 0
        for i in range(len(copy)):
            if copy[i] != [now_one,now_two]:
                if (i== len(copy)-1):
                    copy.append([now_one,now_two])
                    queue.put(copy)
                    break
            else:
                break

    def away(self,now_one,now_two,num,past_list): #水を捨てる
        copy = past_list.copy()
        if num==1:
            now_one = 0
        elif num==2:
            now_two = 0
        for i in range(len(copy)):
            if copy[i] != [now_one,now_two]:
                if (i== len(copy)-1):
                    copy.append([now_one,now_two])
                    queue.put(copy)
                    break
            else:
                break

past_op = []
queue = queue.Queue()
num_one = int(input("Bottle1? ")) 
num_two = int(input("Bottle2? "))
num_Want = int(input("How much do you want? "))

bottle = Bottle(num_one,num_two)
bottle.one = bottle.one_max
past_op.append([bottle.one,bottle.two])
queue.put(past_op)

while (True):
    get_queue = queue.get()
    now = get_queue[-1]
    bottle.one,bottle.two = now[0],now[1]
    if (bottle.one == num_Want or bottle.two == num_Want):
        for i in get_queue:
            print(i)
        break
    if bottle.one > 0 and num_one > bottle.one:
        bottle.put(bottle.one,bottle.two,1,get_queue)
        bottle.pour(bottle.one,bottle.two,1,get_queue)  #oneからtwo
        bottle.away(bottle.one,bottle.two,1,get_queue) #oneを捨てる
    if bottle.one == 0:
        bottle.put(bottle.one,bottle.two,1,get_queue)
    if bottle.one == num_one:
        bottle.pour(bottle.one,bottle.two,1,get_queue)  #oneからtwo
        bottle.away(bottle.one,bottle.two,1,get_queue)

    if bottle.two > 0 and num_two > bottle.two:
        bottle.put(bottle.one,bottle.two,2,get_queue)
        bottle.pour(bottle.one,bottle.two,2,get_queue)  #oneからtwo
        bottle.away(bottle.one,bottle.two,2,get_queue) #oneを捨てる
    if bottle.two == 0:
        bottle.put(bottle.one,bottle.two,2,get_queue)
    if bottle.two == num_two:
        bottle.pour(bottle.one,bottle.two,2,get_queue)  #oneからtwo
        bottle.away(bottle.one,bottle.two,2,get_queue)
    if(queue.empty()):
        print("無理です")
        break

はい、めちゃくちゃ長いです。頑張ればもうちょい短くなると思います。
やってることはBottleクラスでボトル１と２の入っている水量とサイズを定義し、水を入れる関数と水を移す関数、水を捨てる関数を定義します。
あとはキューに今までの水の推移を表したリストを入れていき、ボトル１か２が目的の水量になったらその推移を表したリストを表示して終了です。
もし不可能な組み合わせ(例:6Lと3Lの容器で4Lを測るなど)をしたらいずれキューが空になるのでそうなったら無理と表示して終わりです。

Javascript(+HTML+CSS)のコード

water_jug.html

<!DOCTYPE html>
<html lang="ja">
<head>
    <meta charset="utf-8">
    <title>水量を測る</title>
    <h1>水を測ってみよう！</h1>
    <style type="text/css">
        body{
            text-align: center;
        }
        #table{
            position: relative;
        }
        #table table{
            position:absolute;
            left: 50%;
        }
    </style>
    <script>
        var queue = [];
        var past_op = [];
        function Bottle(one,two){
            this.one = 0;
            this.two = 0;
            this.one_max = one;
            this.two_max = two;
        }
        function show(list_water){
                var table = document.getElementById("table");
                var string = "<table><tr><th>ボトル１</th><th>ボトル２</th></tr><br>";
                for(var i = 0; i<= list_water.length-1; i++){
                    string += "<tr><td>"+list_water[i][0]+"</td><td>"+list_water[i][1]+"</td></tr><br>";
                }
                string +="</table><br>";
                table.innertHTML = string;
        }
        function calc_water(myForm) {

            var bottle = new Bottle(Number(myform.Bottle1.value),Number(myform.Bottle2.value));
            bottle.one = bottle.one_max;
            past_op.push([bottle.one,bottle.two])
            queue.push(past_op)

            while (1){
                var get_queue = queue.shift()
                var now = get_queue[get_queue.length -1]

                bottle.one = now[0];
                bottle.two = now[1];
                if (bottle.one == Number(myform.Want_Water.value) || bottle.two == Number(myform.Want_Water.value)){
                    for (var i=0; i < get_queue.length; i++){
                        console.log(get_queue[i]);
                    }
                    show(get_queue);
                    break;
                }
                if (bottle.one > 0 && bottle.one_max > bottle.one){
                    bottle.put(bottle.one,bottle.two,1,get_queue);
                    bottle.pour(bottle.one,bottle.two,1,get_queue);  //oneからtwo
                    bottle.away(bottle.one,bottle.two,1,get_queue); //oneを捨てる
                }
                if (bottle.one == 0){
                    bottle.put(bottle.one,bottle.two,1,get_queue);
                }
                if (bottle.one == bottle.one_max){
                    bottle.pour(bottle.one,bottle.two,1,get_queue);  //oneからtwo
                    bottle.away(bottle.one,bottle.two,1,get_queue);
                }
                if (bottle.two > 0 && bottle.two_max > bottle.two){
                    bottle.put(bottle.one,bottle.two,2,get_queue);
                    bottle.pour(bottle.one,bottle.two,2,get_queue);  //oneからtwo
                    bottle.away(bottle.one,bottle.two,2,get_queue); //oneを捨てる
                }
                if (bottle.two == 0){
                    bottle.put(bottle.one,bottle.two,2,get_queue);
                }
                if (bottle.two == bottle.two_max){
                    bottle.pour(bottle.one,bottle.two,2,get_queue);  //oneからtwo
                    bottle.away(bottle.one,bottle.two,2,get_queue);
                }
                if (queue.length == 0){
                    cant();
                    break;
                }
            }
        }
        Bottle.prototype.put = function(now_one,now_two,num,past_list){
            var copy = Array.from(past_list);
            if (num ==1){
                now_one = this.one_max;
            }
            else if (num ==2){
                now_two = this.two_max;
            }
            for (var i=0;i< copy.length;i++){
                if (past_list[i].toString() != [now_one,now_two].toString()){
                    if (i==copy.length-1){
                        copy.push([now_one,now_two]);
                        queue.push(copy);
                        break;
                    }
                }
                else{
                    break;
                }
            }
        }
        Bottle.prototype.pour = function(now_one,now_two,num,past_list){
            var copy = Array.from(past_list)
            if (num==1) {
                now_two += now_one;
                if (now_two >= this.two_max){ //あふれる場合,全部は入れない場合
                    now_one = now_two - this.two_max;
                    now_two = this.two_max;
                }
                else {
                    now_one = 0;
                }
            }
            else if (num ==2){ //twoからone
                now_one += now_two;
                if (now_one >= this.one_max){//あふれる場合,全部は入れない場合
                    now_two = now_one - this.one_max;
                    now_one = this.one_max;
                }
                else{
                    now_two = 0;
                }
            }
            for (var i=0;i< copy.length;i++){
                if (past_list[i].toString() != [now_one,now_two].toString()){
                    if (i==copy.length-1){
                        copy.push([now_one,now_two]);
                        queue.push(copy);
                        break;
                    }
                }
                else{
                    break;
                }
            }
        }
        Bottle.prototype.away = function(now_one,now_two,num,past_list){ //水を捨てる
            var copy = Array.from(past_list)
            if (num==1){
                now_one = 0;
            }
            else if (num==2){
                now_two = 0;
            }
            for (var i=0;i< copy.length;i++){
                if (past_list[i].toString() != [now_one,now_two].toString()){
                    if (i==copy.length-1){
                        copy.push([now_one,now_two]);
                        queue.push(copy);
                        break;
                    }
                }
                else{
                    break;
                }
            }
        }  
    </script>
</head>
<body>
    <form id="myform" name="myform"><p>
            Bottle1の大きさ:<input type="text" name="Bottle1"><br>
            Bottle2の大きさ:<input type="text" name="Bottle2"><br>
            欲しい水の量:<input type="text" name="Want_Water"><br> 
            <input type="button" value="計算" 
                name="calc" onclick="calc_water(this.form)"><br>  
    </p></form>
    <div id = "table">
    </div>
</body>
</html>

javascriptは勉強して１週間でクソザコナメクジなのでもっといいコードが書けるかもしれません
やってることはpythonとほとんど同じです。解が見つかったら表を作成するHTMLコードを作成し表示するって感じです。
質問またはこうした方がいいってのがあったら気軽にコメントしてください

Motive

単にスプレッドシートにデータを入力したい場合は gspread ライブラリを使えば sheet.update_acell('A1', "AAA") でよいし、 cellに線や色を付ける、文字のタイプ、文字フォント、または文字寄せなどを設定する場合は gspread-formatting を使って

import gspread_formatting as gsformat

fmt = gsformat.cellFormat(
    backgroundColor=gsformat.color(0.0, 0.5, 0.0),
    textFormat=gsformat.textFormat(foregroundColor=gsformat.color(0, 0, 0)),
    horizontalAlignment='LEFT',
)

gsformat.format_cell_range(sheet, 'A1:A10', fmt)

とすると反映は可能です。

ですが、
A1:A5までcellを結合したい！
思った場合、gspread、gspread-formattingにそれっぽいメソッドはないです。

レイアウトを整えるのは重要事項なのでもっとコードでなんとかしたい！

と思った時は、下記の例はresponseで送るjson形式に似たものを用いる方法をまとめたので記述します。

Method

対処法としてですが gspread の batch_update(dict) メソッドを用います。

code

value.json

{
    "requests" : [
        {
            "mergeCells": {
                "range": {
                    "startRowIndex": 0,
                    "endRowIndex": 1,
                    "startColumnIndex": 0,
                    "endColumnIndex": 5
                },
                "mergeType": "MERGE_ALL"
            }
        }
    ]
}

import configparser
import json
import os
import json
from oauth2client.service_account import ServiceAccountCredentials
import gspread

def main():

    config = configparser.ConfigParser()
    ini_file = os.path.join("./", 'setting.ini')
    config.read(ini_file)

    scope = ["https://spreadsheets.google.com/feeds",
            "https://www.googleapis.com/auth/drive"]
    book_id = config.get("googleSpreadSheet", "book_id")
    path = os.path.join("./", config.get("googleSpreadSheet", "keyfile_name"))
    credentials = ServiceAccountCredentials.from_json_keyfile_name(path, scope)        
    client = gspread.authorize(credentials)
    gfile = client.open_by_key(book_id)

    dst = {}
    with open("value.json", encoding='utf-8') as fin:
        dst = json.load(fin)

    gfile.batch_update(dst)

if __name__ == "__main__":
    main()

これで A1:A5 のセルが結合されました。

Future

batch_update メソッドはcell結合以外にもcopy&pasteやsortもできるみたいで、内部ではPOST https://sheets.googleapis.com/v4/spreadsheets/{spreadsheetId}:batchUpdateにリクエストをかけています。ヒントは(https://developers.google.com/sheets/api/reference/rest/v4/spreadsheets/request)
にあります。

Reference

gspread
gspread-formatting
PythonでGoogleSpreadSheetのセルにコメントを付加したり色や枠線をつける
Basic Formatting
Method: spreadsheets.batchUpdate

　最近ビットコインがまた上昇しているという話ですね。五月あたりに一気に上昇したとかで私のTwitterのTLにいるビットコインの取引している人はめちゃめちゃ喜んでました。そんなわけで今回はブロックチェーンがどのような仕組みになっているかGoとpython(Goとpythonで同じ内容を実装します。所属学部の関係上、pythonが読める人が友達に多いので)で実装しながら説明していきます。最後にまあ最低限ブロックチェーンのシミュレーションといっても差し支えないだろうというもの（主観です）をGoで実装します。

ブロックチェーンとは何か

ブロックチェーン-Wikipedia

ブロックチェーンとは、分散型台帳技術、または、分散型ネットワークである。ビットコインの中核技術（サトシ・ナカモトが開発）を原型とするデータベースである。

とあるように新しい種類のデータベースです。（これは私の考えです。でもたぶんあってます。）よく知られている応用例としてbitcoinやethereumなんかがあります。基本的にbitcoinはお金の記録を保存するもの、ethereumはそれに加えてプログラムの保存をするもののように考えています。ただ、近年はbitcoinにもethereumと同じような機能が実装されてきているらしいです。ブロックチェーンと聞くと「インターネット上にあるお金」を想像する方もいらっしゃると思うんですけど、ブロックチェーンの技術自体はお金とは無関係です。しかしその防御力の高さからかお金の保存や取引などに利用されます。(嘘です。みんながお金をかけても掘りたくなるような実装が組み込まれてるんですけど、それは今回はやりません。)

ブロックチェーンの構造

ハッシュ

ブロックチェーンを理解するのに欠かせない知識の一つにハッシュというものがあります。

ハッシュ-Wikipedia

ハッシュ、ハッシュ値 - データから算出した小さな値。各データを区別・表現する目的に用いる。

>>certutil -hashfile hogehoge.txt SHA256

windowsのコマンドプロンプトでこのように入力するとファイルのハッシュ値を計算することができます。ハッシュ値を算出するための関数をハッシュ関数といいます。ハッシュ関数は出力された値から入力された値を逆算することができないという特徴を持ちます。二つの異なるデータのハッシュ値は互に異なるように設計されているため（必ずしも違う値が出力されるわけではない）、データのハッシュ値が違う場合、互いのデータが違う値であるということが言えます。パスワードの保管などセキュリティが重視される場合、ハッシュ関数にはハッシュの衝突（別のデータのハッシュ値を計算したときに同じ値が出力される現象）ができるだけ少ないSHA(secure hash algorithmの略)と呼ばれる種類のハッシュ関数が使われます。

＊ハッシュ値は5次関数に値を代入したらすぐに答えが求まるのに対して5次方程式は解けないというのとイメージは同じです。（たぶん）

基本構造

ブロックチェーンの基本構造は上の図のように（わかりにくかったらごめんなさい）あるブロックはその前のブロックのハッシュ値をもっているというリンクリストのような構造になっています。先頭は左のブロックです。ブロックチェーンの名前の由来はここです。これで何がうれしいかを説明します。まず、攻撃者が先頭のブロックのデータを書き換えたとします。そうすると次のブロックの中にある先頭のブロックのハッシュ値と、書き換えられた先頭のブロックのハッシュ値が一致しなくなります。そこで攻撃者は矛盾点を消すため二番目のブロックのハッシュ値を書き換えます。すると今度は二番目と三番目のブロックの間に矛盾が生じます。矛盾を消すために攻撃者が三番目のブロックを書き換えます。このように情報を書き換えようと攻撃した場合、書き換えたブロックから一番最後のブロックまですべてのブロックの情報を書き換えなければならないことになります。とりあえずここまでで攻撃するのはいっぱい書き換えなきゃいけないしなんか大変そうだなあってことがわかっていただければよいです。

Goでプログラムを書いてみる

ここまでの説明をプログラムにしてみます。

まずはブロックの構造体とブロックを格納する配列(ブロックチェーン)を定義します。

//後で使うパッケージをimportしとく
package main

import
(
    "fmt"
    "bytes"
    "crypto/sha256"
    "encoding/hex"
    "encoding/json"
)

//ブロックの構造体
type block struct {
    Prevhash [32]byte
}
//都合がいいのでハッシュ値を入れるとブロックを作ってくれる関数も作っとく
//ブロックを作成
func makeBlock(hash [32]byte) block{
    var a block
    a.Prevhash = hash
    return a
}

//ブロックチェーン
//最初のブロックのハッシュ値は0にしとく
var chain = []block{makeBlock([32]byte{})}

次にブロックをjson文字列にする関数と[]byteを渡すとSHA256のハッシュ値を計算してくれる関数を定義

//ブロックをjson文字列に変換する
func convjson(b block) string{
    var buf bytes.Buffer
    byt, _ := json.Marshal(b)
    buf.Write(byt)
    return buf.String()
}

//ハッシュ値を計算する関数
func hash(str []byte)  [32]byte {
    return [32]byte(sha256.Sum256(str))
}

最後にブロックをつなげてく関数をつくって終了！

//マイニング
func mine() string{
    lastBlock := chain[len(chain) - 1]
    lastBlock_json := convjson(lastBlock)
    n := makeBlock(hash([]byte(lastBlock_json)))
    chain = append(chain,n)
    return convjson(n)
}

//ハッシュ値はbyteの配列で保存されています。つまり1byteごとの数字の羅列
//かっちょよく16進数に変換したい場合はこの関数も定義しときましょう。
//[32]byteを[]byteに変換して16進数文字列に変換する
func convhex(b [32]byte) string{
    by := make([]byte,32)
    for i:=0;i<32;i++{
        by[i] = b[i]
    }
    return hex.EncodeToString(by)
}

ホリホリ

func main(){
    //マイニングの回数
    n := 5
    for i:=0;i<n;i++{
        mine()
    }
    fmt.Println(chain)
    /*
    さっきのconvhex関数を定義してたらこっちの方がカッコイイ
    for i:=0;i<len(chain);i++{
        fmt.Println(convhex(chain[i].Prevhash))
    }
     */
}

pythonでプログラムを書いてみる

pythonはGoよりもいくらか簡単にプログラムを書くことができます。すごく短いので説明はコメントアウトしてちょちょっと書きます。さっきと違うこととしては、ブロックを構造体ではなくdictionary(辞書)で作りました。

from hashlib import sha256 # ハッシュ関数
import json #dictionaryをjsonの文字列に変換するため

class blockChain():
    def __init__(self):
        self.chain = [{"prevhash":'0'*64}] # 初期値は0に
    #ディクショナリ、辞書、連想配列（いろんな呼び方があります。）
    #それでブロックを作る関数
    def makeBlock(self,hash):
        b = {} #新しいブロックをdictionaryで
        b["prevhash"] = hash
        return b
    #ブロックののjsonを取得しハッシュ値を計算する
    def hash(self,map):
        data = json.dumps(map) #dictionaryをjsonに変換
        #文字列をエンコードしてハッシュ値を計算して十六進数の文字列に変換
        hash = sha256(data.encode()).hexdigest()
        return hash
    #ブロックをつなげていく
    def mine(self):
        lastBlock = self.chain[-1] # 一番最後のブロックを取得
        hash = self.hash(lastBlock) # ハッシュ値を計算
        nextBlock = self.makeBlock(hash) # 新しいブロックを作る
        self.chain.append(nextBlock) # ブロックをつなげる

c = blockChain()
for i in range(5):
    c.mine()
print(c.chain)

こんな感じで書けます。

ブロックチェーンをより安全にする

上では、ただただ前のブロックのハッシュ値を計算するだけでブロックをつなげていきました。しかし、ハッシュ値を計算するのは一瞬です。これだけの速度でチェーンを伸ばしていけると全て書き換えるのはかんたんです。なので、実際のブロックチェーンは一工夫されています。

上の部分はハッシュ値の計算からブロックの生成までのフローチャートです。前回の実装ではハッシュ値を計算したらそのままブロックを作成しています。しかし、工夫した実装ではハッシュ値に制限をもうけています。性能の良いハッシュ関数は出現する値の確率が一様です。ですので、制限を設けてあげるとそのハッシュ値が出るまで何回も何回も計算しなければいけません。ただし、同じデータについてハッシュ値を計算しても同じ値が出力されるだけです。ですので、ブロックに自由に変更できる要素を追加してその要素をごちゃごちゃ変えて条件を満たすハッシュ値を探していきます。

チェーンは上の図のようになります。これも先頭は左のブロックです。ハッシュ値は入力される値が少し変わるだけで大きく変動するのでごちゃごちゃ値をこねくり回して条件にあうものをさがします。ハッシュ値の制限を厳しくすればするほど合う値を探すのに時間がかかります。こうすることでブロックを改変するのに非常に大きなコストがかかって攻撃が受けにくくなります。

＊ビットコインでは１ブロックつなげるのに平均10分程度かかるようにハッシュ値の制限がかかっています。また、時期によって多くの人が掘る場合やあまり掘られない場合があるので制限がちょうどいいくらいになるように2週間程度に一回制限が修正されます。

Goでプログラムを書いてみる

【変更】まず、構造体にごちゃごちゃ値を変えてハッシュ値を探す値を追加します。この値は一般的にNonceと呼ばれるので変数名にそれを使います。ブロックチェーンの初期値にnonceを追加しときます。

//ブロックチェーン
var chain = []block{makeBlock([32]byte{},0)}

//ブロックチェーンの中身の構造体
type block struct {
    Prevhash [32]byte
    Nonce int //ごちゃごちゃ変える値
}

//ブロックを作成
func makeBlock(hash [32]byte,nonce int) block{
    var a block
    a.Prevhash = hash
    a.Nonce = nonce
    return a
}

【追加】閾値のバイト配列と二つのバイト配列を比較する関数を設定します。aに算出したハッシュ値、bに閾値を入れるとa<bのとき1が返されます。

//256bitのバイト列の比較 a < bのとき1
func comper(a [32]byte,b []byte) int{
    check := make([]byte,32)
    for i:=0;i<32;i++{
        check[i] = a[i]
    }
    return bytes.Compare(b,check)
}

閾値のバイト配列を生成します。指数の計算を今回は32bitで収まる範囲の整数の計算しかしないのでこんな感じで書いてますが、たぶんmathとかのパッケージにあると思うので使うといいと思います。

//指数
func pow(n int,m int) int{
    if m == 1{
        return n
    }else if m == 0 {
        return 1
    }else if m % 2 == 0{
        k := pow(n,m/2)
        return  k*k
    }else {
        return n*pow(n,m-1)
    }
}
//2^nの256bitのバイト列を生成
func makebyte(n int) []byte{
    ans := make([]byte,32)
    //256以上の値が引数に渡されたら制限がないといってよいので
    //とりあえず一番大きい桁に255入れときます。
    //(凝る人はすべてに255入れるループ書いても良いよ)
    if n >= 256{
        ans[0] = byte(255)
        return ans
    }
    amari := n % 8
    syo := n / 8
    check := pow(2,amari)
    ans[31-syo] = byte(check)
    return ans
}

【変更】mine関数のなかで、望んだ範囲のハッシュ値が算出されるまで何度もnonceを変更しながら計算するようにコードを変更します。

//マイニング diffが小さいほど制約が大きい
func mine(diff int) string{
    lastBlock := chain[len(chain) - 1]
    lastBlock_json := convjson(lastBlock)
    lastBlock_hash := hash([]byte(lastBlock_json))
    counter := 0
    limit := makebyte(diff)
    var n block
    for {
        n = makeBlock(lastBlock_hash,counter)
        nowBlock_json := convjson(n)
        nowBlock_hash := hash([]byte(nowBlock_json))
        if comper(nowBlock_hash,limit) == 1{
            break
        }
        counter ++
    }
    chain = append(chain,n)
    return convjson(n)
}

動かしてみます。

func main(){
    //マイニングの回数
    n := 5
    for i:=0;i<n;i++{
        block := mine(240) //難易度を240にしてホリホリ
        fmt.Println(string(block))
    }
}

出力
{"Prevhash":[205,238,133,119,21,139,235,115,177,165,36,115,169,188,168,141,124,99,113,3,149,95,2,232,207,114,60,12,88,69,28,28],"Nonce":280}
{"Prevhash":[0,0,6,23,177,48,240,93,239,0,196,205,228,167,78,52,50,128,27,144,125,8,172,233,27,33,169,127,40,252,55,220],"Nonce":18792}
{"Prevhash":[0,0,93,198,62,162,247,57,110,183,62,209,236,28,195,179,116,88,163,84,1,80,249,49,38,43,253,165,176,14,1,124],"Nonce":7932}
{"Prevhash":[0,0,232,216,132,117,74,174,19,142,116,9,172,50,160,160,64,209,74,47,215,0,182,22,132,178,3,178,166,254,253,70],"Nonce":62563}
{"Prevhash":[0,0,57,142,21,48,208,108,153,110,201,129,204,153,207,154,101,222,46,77,167,82,189,90,50,115,19,20,226,196,53,90],"Nonce":95247}

難易度をかえて掘るのに時間がかかることを確認してみてください。

pythonでプログラムを書いてみる

【変更】さっきのmakeBlock関数にnonce値を加えます。nonce値はごちゃごちゃ変えて欲しいハッシュ値を探す値です。

def makeBlock(self,hash,nonce):
    b = {}
    b["prevhash"] = hash
    b["nonce"] = nonce
    return b

【変更】mine関数の引数にdiffという引数を追加します。これは、ハッシュ値の閾値です。ハッシュ値が2^diffよりも小さい値だった場合はブロックをつなげます。違う場合はnonceを変えてハッシュ値を計算しなおします。

def mine(self,diff):
    lastBlock = self.chain[-1]
    hash = self.hash(lastBlock)
    counter = 0
    while True:
        nextBlock = self.makeBlock(hash,counter)
        if int(self.hash(nextBlock),16) < 2**diff:break
        counter += 1
    self.chain.append(nextBlock)
    return nextBlock

実行！

c = blockChain()
n = 5
for i in range(n):
    block = c.mine(240)
    print(block)

出力
{'prevhash': 'fb197130829d450616b62e431f55fb639aa7ba5989dc82282aeb7b3435856677', 'nonce': 76517}
{'prevhash': '0000c25cb0503255c52ad386cde40b418ec7e0b876dde7cb40039195621c584d', 'nonce': 60994}
{'prevhash': '0000c4180e4cdc8bfaf24ec8738af39800bd08ef9926ab47a56a203194408276', 'nonce': 60232}
{'prevhash': '000030a575effc930f2c4d8529b39ad0e3052853e01b9d3e277dce3cd560877f', 'nonce': 240649}
{'prevhash': '000018c5d3e7bab0564d1b8dce874ed1efd0722c79427709a7ea304537b6cae0', 'nonce': 102298}

ブロックチェーンが安全である理由

ここまでで、新しいブロックを生成するには非常に大きなコストがかかることがわかりました。しかし、これだけでは時間をかけて計算すればすべてのブロックを改変してしまうことができてしまうと思うと思います。ここで、Aさん、Bさん、Cさん、Dさんの四人がSumiSumicoinという暗号通貨をマイニングするときを考えます。四人が掘り始めていくらか時間が経ちました。その時に四人それぞれが掘り進めたブロックのチェーンが下のようになっているとします。

このとき、全員のブロックが違ってはいけないので全員のブロックを統一します。結論から言いますと、統一するために一番長くチェーンをつなげた人のブロックを採用します。今回の場合でいくとBさんのブロックが採用されて、Aさん、Cさん、DさんはBさんのつなげたブロックをコピーします。また、ブロックのつながりに矛盾が生じるもの（ハッシュ値が一致しないもの）は採用しないようにします。このようにすることで、攻撃者がブロックの改ざんをする場合にはほかのまじめにマイニングしている人全員よりも速い速度でブロックをつなげる必要が出てきます。つまり、攻撃者はすべてのマイニングをする人の総和よりも多くの計算リソースを用意する必要が出てくるためそれはほぼ不可能です。私はそれだけ多くのコンピューターを用意して攻撃して採算が取れるとは思いません。まともな人間ならまずやらないでしょう。（個人的な意見です）

この後は複数人でマイニングをした場合を一つのパソコンでシミレーションをするということを行います。ですが、なかなかに長いコードになってしまい、型のないpythonで書くのは辛くなってしまったのでGoのみの実装となります。ご了承ください。

Goで実装してみる

上にも書いたように複数人でマイニングをする場合を実装します。並列処理を実装することであたかも複数人でマイニングをしているように見せかけることができます。

とりあえずこんだけのパッケージをimportしとくと良いです。

import
(
    "math/rand"
    "fmt"
    "reflect"
    "bytes"
    "crypto/sha256"
    "encoding/hex"
    "encoding/json"
    "time"
    "os"
)

まずはブロックの構造を変更します。複数人で掘るということで、マイナーの名前をブロックに記録することにします。また、他人のチェーンが整合性が取れているかどうか判別するためにハッシュ値の制限の閾値を入れてあげます。その変更に合わせてブロックを作成する関数も変更します。

//ブロックチェーンの中身の構造体
type block struct {
    Prevhash [32]byte
    Nonce int
    Difficulty int //難易度
    Name string //マイナーの名前
}

//ブロックを作成
func makeBlock(hash [32]byte,nonce int,difficulty int,name string) block{
    var a block
    a.Prevhash = hash
    a.Nonce = nonce
    a.Difficulty = difficulty
    a.Name = name
    return a
}

今回はブロックチェーンを複数人で掘ります。ですので、人数分のチェーンが必要になります。これに合わせてチェーンを初期化する関数も用意します。

//ブロックチェーン
//var chain = []block{makeBlock([32]byte{},0)}
//チェーンのリスト
var list = make([][]block,4)
//初期化
func initList(){
    for i:=0;i<4;i++{
        list[i] = []block{makeBlock([32]byte{},0,256,"")}
    }
}

チェーンの整合性が取れているかどうかを判別します。これは、他人が悪意のあるチェーンを作った時に自分のチェーンと間違えて置き換えないようにするためです。

//チェーンの整合性をチェック
func check(ls []block) bool{
    //前後のブロックとハッシュ値の関係に矛盾がないかどうか
    for i:=0;i<len(ls) - 1;i++{
        //前のブロックのハッシュ値
        before := hash([]byte(convjson(ls[i])))
        //次のブロックに書いてある前のブロックのハッシュ値
        after := ls[i+1].Prevhash
        if ! reflect.DeepEqual(after,before){
            return false
        }
    }
    //ハッシュ値が2^difficultyよりも小さいかどうか
    for i:=0;i<len(ls);i++{
        //制限のバイト配列
        limit := makebyte(ls[i].Difficulty)
        //現在のブロックのハッシュ値
        nowBlock_hash := hash([]byte(convjson(ls[i])))
        if comper(nowBlock_hash,limit) != 1{
            fmt.Println("yei")
            return false
        }
    }
    return true
}

一番長くて整合性のあるチェーンを返す関数です。ポインタだけ渡されるみたいなことがあると怖いので、最後にコピー走らせてます。

//全てのチェーンを調べて一番長くて整合性のあるチェーンを返す
func consensys() []block{
    //チェーンの長さを入れる。整合性のないチェーンの場合は-1
    chainLength := make([]int,len(list))
    for i:=0;i<len(list);i++{
        if check(list[i]){
            chainLength[i] = len(list[i])
        } else {
            chainLength[i] = -1
        }
    }
    maxLength := 0
    maxIndex := 0
    for i:=0;i<len(chainLength);i++{
        if chainLength[i] > maxLength{
            maxLength = chainLength[i]
            maxIndex = i
        }
    }
    ansChain := make([]block,maxLength)
    for i:=0;i<maxLength;i++{
        ansChain[i] = list[maxIndex][i]
    }
    return ansChain
}

mine関数も変えます。まず、1から順に探していたnonce値を乱数にします。これは、複数の人間が全く同じようにnonce値を探していたらちょっと早く掘り始めた人の一人勝ちになってしまうからです。（たぶん）実験はしてません。たぶんそうです。たぶん．．．
また、ひとブロック掘るごとに一番長いチェーンを自分のチェーンに置き換えます。この実装だと長いチェーンが自分だとしても自分のチェーンをコピーして置き換えるという操作をしているので凝る人は高速化できると思います。引数のindexはさっき作ったチェーンのリストのインデックスを指定してください。

//マイニング diffが小さいほど制約が大きい indexはチェーンのリストの何番目か
func mine(diff int,index int,name string) string{
    lastBlock := list[index][len(list[index]) - 1]
    lastBlock_json := convjson(lastBlock)
    lastBlock_hash := hash([]byte(lastBlock_json))
    limit := makebyte(diff)
    var n block
    for {
        //                          乱数にする↓
        n = makeBlock(lastBlock_hash,rand.Int(),diff,name)
        nowBlock_json := convjson(n)
        nowBlock_hash := hash([]byte(nowBlock_json))
        if comper(nowBlock_hash,limit) == 1{
            break
        }
    }
    //一番長いチェーンを自分の所に置き換える
    list[index] = append(list[index],n)
    list[index] = consensys()
    return convjson(n)
}

最後にホリホリするためのmain関数を書きます。A、B、C、Dの四人でマイニングをします。それぞれmine関数に名前とインデックスを指定します。また、もう一つ5秒おきに最長のチェーンの状態を出力してくれるものも並列で実行しときます。

func main(){
    initList()
    //Aさん
    go func() {
        for {
            mine(235,0,"A")
        }
    }()
    //Bさん
    go func() {
        for {
            mine(235,1,"B")
        }
    }()
    //Cさん
    go func() {
        for {
            mine(235,2,"C")
        }
    }()
    //Dさん
    go func() {
        for {
            mine(235,3,"D")
        }
    }()

    //五秒おきに最長チェーンの状態を出力
    go func() {
        for {
            time.Sleep(5*time.Second)
            ls := consensys()
            for i:=0;i<len(ls);i++{
                fmt.Println(ls[i])
            }
        }
    }()
    //四人で5分マイニングした後最長チェーンの状態をファイルに書き込む
    time.Sleep(300*time.Second)
    //                     ファイル名を指定↓
    file, err := os.Create(`directory\\hogehoge.txt`)
    if err != nil {
        // Openエラー処理
    }
    defer file.Close()

    ls := consensys()
    for i:=0;i<len(ls);i++{
        output := convjson(ls[i]) + "\r\n"
        file.Write(([]byte)(output))
    }
}

ここまで出来たら実行してあげましょう。五分後には自作ブロックチェーンの結果がファイルに出力されているはずです。私のPC(2コア4スレッド2.50GHz - 2.71GHz)では235がちょうどよい難易度だったのでそれになってまずが、もっといいプロセッサが積んであるコンピューターではもう少し数字を下げてあげてもいいかもしれません。

私が出力したやつ

このソースコード

最後に

ここまで作ってきましたが、やはりブロックチェーンの構造は面白いですね。これでとりあえず変更されにくいチェーンを組めたので、あとはブロックの中にどんな機能を入れるかによってさまざまな性質のブロックチェーンを作ることができます。世の中にたくさんの暗号通貨が存在しているのはそのためです。（ほかにさして知識のない人が自分の通貨を作りたいと言って作る人もいますが...）なんにせよ、ブロックチェーンは未来ある技術なので今後の展開が楽しみです。あと、私はこれからテスト期間に入るので勉強頑張ります。

はじめに

相変わらず競馬関係のプログラムを書いています。
これまではHeroku上で開発していたのですが、仕事でAWSを使っていることもあり、個人開発でもAWSを使い始めました。

Herokuだと毎日決まった時間に動作する処理を書く場合、Heroku Schedulerを使う形になりますが、AWSだとLambdaをCloudWatchイベントでキックする形…とはいうもののひとまとめになっている情報が見当たらなかったので、やってみたことを残しておきたいと思います。

※2019年7月時点の情報です。AWSはサービス改善のスピードが速すぎてWebにまとめた情報があっという間に陳腐化するので、そこはご注意を…

AWS SAMのテンプレート

AWS SAMを使ってみたのですが、元になるテンプレートを作成するためにsam initコマンドを使うとAPIゲートウェイからキックされる処理のテンプレートになります。

定時実行する処理用のテンプレートは…ということで調べたところ、以下にありました。
CloudWatch イベント アプリケーションの AWS SAM テンプレート

だけど、これは定時実行ではなくて定期実行、かつ今回必要ないようなリソースも使っていたりだったので、定期実行する最低限のテンプレートを作りました。

template.yaml

AWSTemplateFormatVersion: '2010-09-09'
Transform: AWS::Serverless-2016-10-31
Description: >
  notify_hitokuchi

  Sample SAM Template for notify_hitokuchi

Globals:
  Function:
    Timeout: 10

Resources:
  NotifyHitokuchiFunction:
    Type: AWS::Serverless::Function
    Properties:
      CodeUri: notify_hitokuchi/
      Handler: app.lambda_handler
      Runtime: python3.7
      Events:
        NotifyHitokuchi:
          Type: Schedule
          Properties:
            Schedule: cron(30 12 * * ? *)

一番下の行のScheduleの行が定時実行する時刻の指定。式の記述方法は以下を参照。
ルールのスケジュール式

ハンドラでもらえるイベント情報の時刻

先程のテンプレートだと、notify_hitokuchi/app.py内のlambda_handler関数が実行時に呼び出される関数になります。

app.py

def lambda_handler(event, context):

eventには辞書形式でイベント情報が入ってきます。
今回はCloudWatchイベントでキックされた時刻が欲しかったのですが、実行時刻はevent['time']にISO 8601形式で入ってきます。
それをそのまま「datetime.fromisoformat()」関数に渡せばdatetimeが入手できる…と思いきや、「2019-07-10T12:30:00Z」という形だと「datetime.fromisoformat()」は処理できないようでして…
以下のようなコードでdatetimeを入手しました。

# event['time']に'2019-07-10T12:30:00Z'形式の文字列が入ってくる
datetime_str = event['time'].replace('Z', '+00:00')

# '2019-07-10T12:30:00+00:00'形式だとfromisoformatで処理できる
datetime_utc = datetime.fromisoformat(date_str)

pip等でインストールしたライブラリ

pip等でインストールしたライブラリをAWS Lambdaで使う場合、ライブラリを1箇所にまとめてZIPで固めて…という話をよく見かけますが、samを使っている分にはsam buildで済むようです。

$ sam build
2019-07-19 07:42:34 Building resource 'NotifyHitokuchiFunction'
2019-07-19 07:42:34 Running PythonPipBuilder:ResolveDependencies
2019-07-19 07:42:44 Running PythonPipBuilder:CopySource

Build Succeeded

Built Artifacts  : .aws-sam/build
Built Template   : .aws-sam/build/template.yaml

Commands you can use next
=========================
[*] Invoke Function: sam local invoke
[*] Package: sam package --s3-bucket <yourbucket>

↓↓↓2019-07-20追記↓↓↓
sam buildの前にrequirements.txtをCodeUriで指定したパスに設置が必要です。
自分はpipenvを使っているので、以下でrequirements.txtを生成しました。

$ $ pipenv lock -r > notify_hitokuchi/requirements.txt 
$ cat notify_hitokuchi/requirements.txt 
-i https://pypi.org/simple
beautifulsoup4==4.7.1
bs4==0.0.1
certifi==2019.6.16
chardet==3.0.4
decorator==4.4.0
html5lib==1.0.1
idna==2.8
mojimoji==0.0.9.post0
py==1.8.0
python-dotenv==0.10.3
requests==2.22.0
retry==0.9.2
six==1.12.0
soupsieve==1.9.2
urllib3==1.25.3
webencodings==0.5.1

また、上記requirements.txtに記載のあるmojimojiのようにC/C++で作られているライブラリは動作環境でコンパイルする必要がありますが、sam buildに--use-containerオプションを付けることでdockerコンテナ上でビルドが行われ、Lambda上で動くライブラリを作ることができます。素晴らしい!

$ sam build --use-container
2019-07-20 09:43:51 Starting Build inside a container
2019-07-20 09:43:51 Building resource 'NotifyHitokuchiFunction'

Fetching lambci/lambda:build-python3.7 Docker container image......
2019-07-20 09:43:55 Mounting /Users/masaminh/develop/notify_hitokuchi_aws/notify_hitokuchi as /tmp/samcli/source:ro,delegated inside runtime container

Build Succeeded

Built Artifacts  : .aws-sam/build
Built Template   : .aws-sam/build/template.yaml

Commands you can use next
=========================
[*] Invoke Function: sam local invoke
[*] Package: sam package --s3-bucket <yourbucket>

Running PythonPipBuilder:ResolveDependencies
Running PythonPipBuilder:CopySource

↑↑↑2019-07-20追記↑↑↑

ローカルでの動作確認

まずハンドラに渡されるイベント情報を作っておく必要がありますが、その元はsam local generate-eventを使って取得できます。
今回はCloudWatchのイベントですが、以下で作成しました。

$ sam local generate-event cloudwatch scheduled-event > event.json
$ cat event.json
{
  "id": "cdc73f9d-aea9-11e3-9d5a-835b769c0d9c",
  "detail-type": "Scheduled Event",
  "source": "aws.events",
  "account": "",
  "time": "1970-01-01T00:00:00Z",
  "region": "us-east-1",
  "resources": [
    "arn:aws:events:us-east-1:123456789012:rule/ExampleRule"
  ],
  "detail": {}
}

実際にローカルで動作確認を行う際には、

$ sam local invoke -e event.json NotifyHitokuchiFunction
2019-07-19 08:03:10 Invoking app.lambda_handler (python3.7)
2019-07-19 08:03:10 Found credentials in shared credentials file: ~/.aws/credentials

Fetching lambci/lambda:python3.7 Docker container image......

という感じで呼び出すことができます。Dockerコンテナ上で動くので起動するのに少し時間がかかります。

デプロイ

まずパッケージを作成。

$ sam package --s3-bucket バケット名 --output-template-file package.yaml

Successfully packaged artifacts and wrote output template to file package.yaml.
Execute the following command to deploy the packaged template
aws cloudformation deploy --template-file /Users/masaminh/develop/notify_hitokuchi_aws/package.yaml --stack-name <YOUR STACK NAME>

引き続き、メッセージに従ってデプロイ...

$ aws cloudformation deploy --template-file package.yaml --stack-name notify-hitokuchi

Waiting for changeset to be created..

Failed to create the changeset: Waiter ChangeSetCreateComplete failed: Waiter encountered a terminal failure state Status: FAILED. Reason: Requires capabilities : [CAPABILITY_IAM]

メッセージのままだと上記のように怒られるので、capabilitiesにCAPABILITY_IAMを設定します。

$ aws cloudformation deploy --template-file package.yaml --stack-name notify-hitokuchi --capabilities CAPABILITY_IAM

Waiting for changeset to be created..
Waiting for stack create/update to complete
Successfully created/updated stack - notify-hitokuchi

デプロイできました。

ログ出しに関して

app.py

import logging

logger = logging.getLogger()
logger.setLevel(logging.INFO)

def lambda_handler(event, context):
    logger.info('start: event.time=%s', event['time'])

こんな感じで普通にloggerから出力すれば、CloudWatchで見られます。

[INFO] 2019-07-18T12:30:28.336Z 51865e95-cacd-4f9d-80e9-6e1f148357b9 start: event.time=2019-07-18T12:30:00Z

おわりに

今回は競馬的なことは書かずにノウハウ的なところのみでした。(Lambda関数名からして競馬系の処理だなぁという感じではありますが・笑)
引き続き競馬系の個人開発はやってますので、またそのうち競馬開発系の記事を書きたいな〜と思います(笑)

はじめに

この記事は古川研究室 Workout_calendar ４日目の記事です。
本記事は古川研究室の学生が学習の一環として書いたものです。

本記事では重回帰分析について述べます。
重回帰分析を扱った時にプログラムは回ったのに予測値がめちゃくちゃ！ってことはよくありますよね。その原因の一つに学習データの型(サイズ)があるかもしれません。重回帰では入力の次元数よりもデータ数が少ないと、出力値に大きな差が出ることがあります。
本記事では重回帰分析のアルゴリズムについて詳しく述べた後、入力が高次元かつデータ数が少ない場合に何が起こっているのか、またその解決法について述べます。

重回帰分析とは

問題設定

観測データとしてD次元のベクトルの入力$\mathbf x\in \mathbb R^D$と出力$y\in \mathbb R$のペア$\left\{(\mathbf{x}_n,y_n)\right\} _{n=1}^N $から$y=f(\bf{x})$となるような関数$f$を推定すること

重回帰分析を解く

$\begin{equation}
f(\mathbf{x})=w^{(1)}x^{(1)}+w^{(2)}x^{(2)}+...+w^{(D)}x^{(D)} \tag{1}
\end{equation}$
とモデルを定義します。
$\begin{equation}
f(\mathbf{x})=\sum_{d=1}^{D}{w^{(d)}}{x^{(d)}} \tag{2}
\end{equation}$
$\mathbf{w}=({w^{(1)}},{w^{(2)}},...,{w^{(D)}})^T,\mathbf{x}=(x^{(1)},x^{(2)},...x^{(D)})^T$とおくと
$\begin{equation}
f(\mathbf{x})=\mathbf{w}^T\mathbf{x} \tag{3}
\end{equation}$

観測データの出力とモデルの出力の二乗誤差の総和$J_\mathrm{LS}(\mathbf{w})$が最小となる$\mathbf{w}$を求めていきます。

\begin{align}
J_\mathrm{LS}(\mathbf{w})&=\sum_{n=1}^{N}(y_n-f(\mathbf{w}))\\
&=\sum_{n=1}^{N}(y_n-\mathbf{w}^T\mathbf{x}_n)\\
&=\sum_{n=1}^{N}(y_n-\mathbf{x}_n^T\mathbf{w})
\end{align}

\mathbf{y}=(y_1,y_2,...,y_N)^T,\mathbf{X}=\begin{pmatrix}x_1^{(1)} & x_1^{(2)} & \cdots & x_1^{(D)} \\
x_2^{(1)} & x_2^{(2)} & \cdots & x_2^{(D)} \\
\vdots & \vdots & \ddots & \vdots \\
x_N^{(1)} & x_N^{(2)} & \cdots & x_N^{(D)}&\end{pmatrix}

とおくと

\begin{eqnarray}
    J_\mathrm{LS}(\mathbf{w})&=&(\mathbf{y}-\mathbf{X}\mathbf{w})^T(\mathbf{y}-\mathbf{X}\mathbf{w})\\
    &=&(\mathbf{y}^T-\mathbf{w}^T\mathbf{X}^T)(\mathbf{y}-\mathbf{X}\mathbf{w})\\
    &=&\mathbf{y}^T\mathbf{y}-\mathbf{y}^T\mathbf{X}\mathbf{w}-\mathbf{w}^T\mathbf{X}^T\mathbf{y}+\mathbf{w}^T\mathbf{X}^T\mathbf{X}\mathbf{w}\\
    &=&\mathbf{y}^T\mathbf{y}-2\mathbf{w}^T\mathbf{X}^T\mathbf{y}+\mathbf{w}^T\mathbf{X}^T\mathbf{X}\mathbf{w}\\
\end{eqnarray}

上式を$\mathbf{w}$について偏微分し0とおくと

\begin{eqnarray}
    \frac{\partial}{\partial \mathbf{w}}J_\mathrm{LS}&=&\frac{\partial}{\partial \mathbf{w}}(\mathbf{y}^T\mathbf{y}-2\mathbf{w}^T\mathbf{X}^T\mathbf{y}+\mathbf{w}^T\mathbf{X}^T\mathbf{X}\mathbf{w})\\
    0&=&0-2\mathbf{X}^T\mathbf{y}+\{\mathbf{X}^T\mathbf{X}+(\mathbf{X}^T\mathbf{X})^T\}\mathbf{w}\\
    2\mathbf{X}^T\mathbf{X}\mathbf{w}&=&2\mathbf{X}^T\mathbf{y}\\
\mathbf{X}^T\mathbf{X}\mathbf{w}&=&\mathbf{X}^T\mathbf{y}\\
\mathbf{w}&=&(\mathbf{X}^T\mathbf{X})^{-1}\mathbf{X}^T\mathbf{y}
\end{eqnarray}

となりパラメータ$\mathbf{w}$を求めることができました！

具体的な例を当てはめてみる

アイスクリーム屋の一日の売上をその日の気温と湿度より予測する重回帰問題を考えてみます。

# 必要なライブラリのimport
import numpy as np
import random
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# 入力データセット
N = 20
a1 = 1
a2 = 15
temp = np.random.uniform(22,35,N)
humi = np.random.uniform(0.0,0.9,N)
sales = a1*temp + a2*humi + np.random.rand(N)*2 - 1

# キャンバス作成
fig = plt.figure(figsize=(9, 9))
ax = fig.add_subplot(1, 1, 1, projection='3d')
ax.set_title("example data: sales of ice")
ax.set_xlabel("humidity")
ax.set_ylabel("temperature")
ax.set_zlabel("sales of ice")
ax.scatter(temp, humi, sales,s=30, c='red')

# データセットを行列に
x_1 =temp [:,None]
x_2 =humi [:,None]
y = sales
X = np.concatenate([x_1,x_2],axis=1)

# 重回帰の解
w = np.linalg.inv(X.T@X)@X.T@y

new_samples_num = 30
new_temp,new_humi=np.meshgrid(np.linspace(22,35,new_samples_num),np.linspace(0.0,0.9))
new_temp=new_temp[:,:,None]
new_humi=new_humi[:,:,None]
new_X=np.concatenate([new_temp, new_humi], axis=2)
new_y=new_X[:,:,0]*w[0]+new_X[:,:,1]*w[1]

fig=plt.figure(1)
ax=Axes3D(fig)
ax.set_xlabel('temperature')
ax.set_ylabel('humidity')
ax.set_zlabel('sales')
ax.set_title('regression by linear model')
ax.scatter(X[:,0],X[:,1],y,c='red',label='trainig data $(x_{n1},x_{n2},y_n)$')
ax.plot_wireframe(new_X[:,:,0], new_X[:,:,1], new_y, color='b',label='estimated function $f_w(x)$')
plt.legend(fontsize=12)
plt.show()

気温と湿度からその日のアイスの売上を予測することができるようになりました！

入力の次元数がデータ数よりも多い場合

以下の重回帰モデルの解を見てみると、重回帰では$\mathbf{X}^T\mathbf{X}$が逆行列を持つ、すなわち正則であることが前提としてあることが分かります。

\begin{eqnarray}
\mathbf{w}&=&(\mathbf{X}^T\mathbf{X})^{-1}\mathbf{X}^T\mathbf{y}
\end{eqnarray}

ここで入力の次元数$D$がデータ数$N$より多い場合を考えてみましょう。$(D>N)$
$\mathbf{X}$を特異値分解します。

\begin{eqnarray}
\mathbf{X}=U\Sigma V^T
\end{eqnarray}

ただし$\Sigma$は$N \times D$の行列で以下に示す性質を持っています。

\begin{eqnarray}
\Sigma &=& [\Delta,\mathbf{0}] \\
\Delta &=& {\rm diag}(\sigma_1, \dots, \sigma_N)
\end{eqnarray}

このとき

\begin{eqnarray}
\mathbf{X}^T\mathbf{X}&=&(U\Sigma V^T)^TU\Sigma V^T\\
&=&V \Sigma^T U^T U \Sigma V^T\\
&=&V \Sigma^T \Sigma V^T\\
&=&V\begin{pmatrix} \Delta ^2 & \mathbf{0} \\
\mathbf{0} & \mathbf{0}
\end{pmatrix}V^T
\end{eqnarray}

となり$\mathbf{X}^T\mathbf{X}$は非正則となります。数学的には、「逆行列は求められない」で終わりなのですが、数値計算上では$\Sigma^T \Sigma$の要素が常識のように完全に$0$になるのではなく、$0$に近い非常に小さな値が入ります。それに対して無理やり逆行列を計算することになるので$\mathbf{w}$の要素にはとても大きな値が入ることになり、そのモデルから導かれた出力間には大きな差がでてしまいます。

解決法

$\mathbf{X}^T\mathbf{X}$が非正則である時にリッジ回帰という手法が使えます。(他にも過学習を防ぐ場合などにも使われます)
リッジ回帰では重みベクトル$\mathbf{w}$の大きさにペナルティを加えます。($\alpha>0$)

\begin{eqnarray}
    J_\mathrm{LS}(\mathbf{w})&=&\|\mathbf{y}-\mathbf{Xw}\|^2+\alpha\|\mathbf{w}\|^2
\end{eqnarray}

同様に$\mathbf{w}$で微分して0とおくと

\begin{eqnarray}
    \mathbf{w}=(\mathbf{X}^T\mathbf{X}+\alpha I)^{-1}\mathbf{X}^T\mathbf{y}
\end{eqnarray}

となりますね。
$\mathbf{X}$の特異値分解を考えると

\begin{eqnarray}
    \mathbf{X}&=&U\Sigma V^T\\
    \mathbf{X}^T\mathbf{X}+\alpha I &=& (U\Sigma V^T)^T U\Sigma V^T+\alpha I\\
    &=&V\Sigma^2V^T+\alpha VIV^T\\
    &=&V(\Sigma^2+\alpha I)V^T
\end{eqnarray}

となります。
$\mathbf{X}^T\mathbf{X}+\alpha I$の特異値は$\mathbf{X}$の特異値の2乗に正の値を足したものとなるため全て正となるので…

\begin{eqnarray}
    \rm{det}(\mathbf{X}^T\mathbf{X}+\alpha I)>0
\end{eqnarray}

となり、$\mathbf{X}^T\mathbf{X}+\alpha I$は正則となって逆行列を持ちます！

まとめ

重回帰分析において入力の次元数に対してデータ数が足りない場合の解決法とそれで解決される理由を示しました。
皆さんも問題が発生した際は自分が扱っている手法が何を前提としているかを考えるのもいいかもしれません。

参考文献：パターン認識と機械学習(上・下)

Pythonを使ってベクトルをL2正規化（normalization）する方法が色々あるのでまとめます。
※L2正則化（regularization）= Ridgeではありません。

L2正規化とは

ベクトル$x = (x_1,x_2,...,x_n)$に対応するL2正規化は以下のように定式化されます。

x_{L2\,norm} = \frac{x}{||x||_2}

ここで、$||x||_2$は、以下の式で求まる$x$のL2ノルムです。

||x||_2 = \sqrt{(\sum_i x_i^2)} = \sqrt{x_1^2 + x_2^2 + ... + x_n^2}

ベクトル$x$をL2正規化すると、長さが1のベクトルになります。

２次元空間で考えた場合、この操作は任意の２次元ベクトルを原点を中心とした半径1の円上へのベクトルに変換していると考えることができます。

３次元の場合は原点を中心とした半径1の球面上へのベクトルに変換します。

４次元以上はイメージが難しいですが、同じように原点を中心とする半径1の超球面上へのベクトルに変換すると考えれば良いと思います。

pythonでの計算方法

上記のようなL2正規化ですが、pythonでは色々と計算方法があるようなので以下にまとめます。
他にもあると思うので、見つけたら追記します。

直接計算する

定義に従って計算すれば良いです。

# 定義に従って計算
x = np.array([1,2,3,4,5])
x_l2_norm = sum(x**2)**0.5
x_l2_normalized = x / x_l2_norm
# 検算
np.linalg.norm(x_l2_normalized,ord=2)
# 1.0

numpy.linalg.norm

numpy.linalg.normを使って計算することも可能です。
こいつはベクトルxのL2ノルムを返すので、L2ノルムを求めた後にxを割ってあげる必要があります。

# numpy.linalg.normを使う
x = np.array([1,2,3,4,5])
x_l2_norm = np.linalg.norm(x,ord=2)
x_l2_normalized = x / x_l2_norm
# 検算
np.linalg.norm(x_l2_normalized,ord=2)
# 1.0

sklearn.preprocessing.normalize

みんな大好きscikit-learnでも計算可能です。
２次元の場合はaxisで軸の指定します。axis=0で列方向、axis=1(デフォルト)で行方向に正規化を行います。

# sklearn.preprocessing.normalizeを使う
from sklearn.preprocessing import normalize
x = [[1,2,3,4,5],]
x_l2_normalized = normalize(x,norm='l2')
x_2d = [[1,2,3,4,5],[6,7,8,9,10]]
x_2d_l2_normalized = normalize(x_2d,norm='l2')
x_2d_l2_normalized2 = normalize(x_2d,norm='l2',axis=0)
# 検算
print(np.linalg.norm(x_l2_normalized,ord=2))
print(np.linalg.norm(x_2d_l2_normalized[0,:],ord=2))
print(np.linalg.norm(x_2d_l2_normalized[1,:],ord=2))
print(np.linalg.norm(x_2d_l2_normalized2[:,0],ord=2))
print(np.linalg.norm(x_2d_l2_normalized2[:,1],ord=2))
print(np.linalg.norm(x_2d_l2_normalized2[:,2],ord=2))
print(np.linalg.norm(x_2d_l2_normalized2[:,3],ord=2))
print(np.linalg.norm(x_2d_l2_normalized2[:,4],ord=2))
# 1.0
# 1.0
# 1.0
# 1.0
# 1.0
# 1.0
# 1.0
# 0.9999999999999999

keras.backend.l2_normalize

kerasにもあります。backendはtensorflowです。
(tf.keras.backend.l2_normalizeというものもありました)

# keras.backend.l2_normalizeを使う
from keras import backend as K
x = [1.,2.,3.,4.,5.]
x_2d = [[1.,2.,3.,4.,5.],[6.,7.,8.,9.,10.]]
x_l2_normalized = K.l2_normalize(x,axis=0)
x_2d_l2_normalized = K.l2_normalize(x_2d,axis=0)
x_2d_l2_normalized2 = K.l2_normalize(x_2d,axis=1)
with tf.Session() as sess:
    sess.run(tf.global_variables_initializer())
    x_l2_normalized = sess.run(x_l2_normalized)
    x_2d_l2_normalized = sess.run(x_2d_l2_normalized)
    x_2d_l2_normalized2 = sess.run(x_2d_l2_normalized2)
# 検算
print(np.linalg.norm(x_l2_normalized,ord=2))
print(np.linalg.norm(x_2d_l2_normalized[:,0],ord=2))
print(np.linalg.norm(x_2d_l2_normalized[:,1],ord=2))
print(np.linalg.norm(x_2d_l2_normalized[:,2],ord=2))
print(np.linalg.norm(x_2d_l2_normalized[:,3],ord=2))
print(np.linalg.norm(x_2d_l2_normalized[:,4],ord=2))
print(np.linalg.norm(x_2d_l2_normalized2[0,:],ord=2))
print(np.linalg.norm(x_2d_l2_normalized2[1,:],ord=2))
#1.0000001
#1.0
#0.99999994
#1.0
#0.99999994
#1.0
#1.0000001
#1.0

tf.math.l2_normalize

※tensorflow 1.11まではtf.nn.l2_normalizeとして提供されていたようです。

import tensorflow as tf

x = tf.Variable([1,2,3,4,5], dtype=tf.float32)
x_2d = tf.Variable([[1,2,3,4,5],[6,7,8,9,10]], dtype=tf.float32)
x_l2_normalized = tf.math.l2_normalize(x,axis=0)
x_2d_l2_normalized = tf.math.l2_normalize(x_2d,axis=0)
x_2d_l2_normalized2 = tf.math.l2_normalize(x_2d,axis=1)

with tf.Session() as sess:
    sess.run(tf.global_variables_initializer())
    x_l2_normalized = sess.run(x_l2_normalized)
    x_2d_l2_normalized = sess.run(x_2d_l2_normalized)
    x_2d_l2_normalized2 = sess.run(x_2d_l2_normalized2)
# 検算
print(np.linalg.norm(x_l2_normalized,ord=2))
print(np.linalg.norm(x_2d_l2_normalized[:,0],ord=2))
print(np.linalg.norm(x_2d_l2_normalized[:,1],ord=2))
print(np.linalg.norm(x_2d_l2_normalized[:,2],ord=2))
print(np.linalg.norm(x_2d_l2_normalized[:,3],ord=2))
print(np.linalg.norm(x_2d_l2_normalized[:,4],ord=2))
print(np.linalg.norm(x_2d_l2_normalized2[0,:],ord=2))
print(np.linalg.norm(x_2d_l2_normalized2[1,:],ord=2))
#1.0000001
#1.0
#0.99999994
#1.0
#0.99999994
#1.0
#1.0000001
#1.0

おわりに

いろいろありますねー。用途に合わせて使い分けていきましょー（？）
（pythonは変数のスコープどうなってんすか。。。）

CNNの畳み込みの計算をpythonで実装してみたらim2colになっていた話。

im2colって？

かの有名な「ゼロから作るdeep learning」の畳み込み演算に出てくる関数。（自分はこの本持ってませんが…）検索するとこんなのが出てきます。

def im2col(input_data, filter_h, filter_w, stride=1, pad=0):
    N, C, H, W = input_data.shape
    out_h = (H + 2*pad - filter_h)//stride + 1
    out_w = (W + 2*pad - filter_w)//stride + 1

    img = np.pad(input_data, [(0,0), (0,0), (pad, pad), (pad, pad)], 'constant')
    col = np.zeros((N, C, filter_h, filter_w, out_h, out_w))

    for y in range(filter_h):
        y_max = y + stride*out_h
        for x in range(filter_w):
            x_max = x + stride*out_w
            col[:, :, y, x, :, :] = img[:, :, y:y_max:stride, x:x_max:stride]

    col = col.transpose(0, 4, 5, 1, 2, 3).reshape(N*out_h*out_w, -1)
    return col

img2colでも引っ掛かりますが、多分同じ内容です。

単純な畳み込み計算

単純にfor文で畳み込みを計算してみます。
入力$(224,224,8)$を出力$(224,224,16)$に変換するサイズ$(3,3)$の二次元畳み込みを考えるとします。
畳み込みをforループの力業で無理に計算しようとすれば以下のような計算をすればよいことになります。

import numpy as np
import itertools

np.random.seed(100)

#(224,224,8)=>(224,224,16)
input = np.ones((224,224,8))
output = np.zeros((224,224,16))
weight = np.random.random_sample((3,3,8,16))

input = np.pad(input, [(1,1),(1,1),(0,0)], 'constant')  # (224,224,8)=>(226,226,8)
x = np.zeros((224,224,8,16))

# case simple conv cal
for ch1, ch2, i, j, k, l in itertools.product(range(8), range(16), range(224), range(224), range(3), range(3)):
    x[i,j,ch1,ch2] += input[i+k,j+l,ch1] * weight[k,l,ch1,ch2]

output = np.sum(x, axis=(2))

output1= [[[17.6167424  16.46092477 16.71833326 ... 15.7476869  16.24364858
   15.14112872]...
time= 74.62933945655823[s]

単純な式で記述出来ますが、この多重for文では計算に非常に時間がかかります。
この単純な畳み込み計算は現実的ではありません。

一次元畳み込みの場合

最初は簡単な一次元畳み込みを代わりに考えてみます。
ここで$*$は畳み込み演算子です。するとこの式は以下のように記述できます。
この整理した$z_i$の要素はZeroPaddingした$(0,x_1,x_2,...,x_{n},0)$からひとつずつずらしながら取り出す長さ$n$の配列になります。
このようにZeroPaddingした配列をずらして得た配列をあらかじめ用意してやり、これに畳み込みfilter要素を行列計算してやると多重for文よりも高速になります。

\begin{align}
y &= A*x \\
y &= (A_1,A_2,A_3) * (x_1, x_2, x_3, ... x_m, ... ,x_n)\\
&=(A_2x_1+A_3x_2, A_1x_1+A_2x_2+A_3x_3, ... \\
&\qquad A_1x_{m-1}+A_2x_{m}+A_3x_{m+1},...\\
&\qquad A_1x_{n-1}+A_2x_{n} )\\
&=A_1(0,x_1,x_2,...,x_{n-2},x_{n-1})\\
&+A_2(x_1,x_2,x_3...,x_{n-1},x_{n})\\
&+A_3(x_2,x_3,x_4...x_{n},0)\\
&=(A_1,A_2,A_3)\cdot(z_1,z_2,z_3)
\end{align}

二次元の畳み込みの場合もこれと同様に考える事ができます。

実装

以下はZeroPaddingしたinputをずらしてxの配列を取り、それを変形して重みとの行列積計算を行います。ここでfilter要素数9×入力チャンネル数8=内積計算の長さ$72$になってます。
念のため範囲を全部書けば$(224,224)$をZeroPaddingした$(226,226)$の配列の内$(0:224,0:224),$$(1:225,0:224),$$(2:226,0:224),$$(0:224,1:225),$$(1:225,1:225),$$(2:226,1:225),$$(0:224,2:226),$$(1:225,2:226)$$,(2:226,2:226)$の9通りの範囲になります。

input = np.pad(input, [(1,1),(1,1),(0,0)], 'constant')  # (224,224,8)=>(226,226,8)

# case im2col
x = np.zeros((224,224,8,3,3))
weight = weight.transpose((2,0,1,3))               # (3,3,8,16) => (8,3,3,16)
weight = np.reshape(weight,(72,16))                # (8,3,3,16) => (72,16)

for i, j in itertools.product(range(3), range(3)):
    x[:,:,:,i,j] = input[i:224+i,j:224+j]

x = np.reshape(x,(224*224,72))                     # (224,224,8,3,3) => (224*224,72)
x = np.dot(x,weight)                               # (224*224,72)*(72,16) => (224*224,16)
output = np.reshape(x,(224,224,16))                # (224*224,16)    => (224,224,16)

output2= [[[17.6167424  16.46092477 16.71833326 ... 15.7476869  16.24364858
   15.14112872]...
time= 0.0625007152557373[s]

この場合、計算時間はかなり短くなりました。また、計算結果も単純な畳み込み計算結果と合ってます。
ここでこの内容を探してみると前述のim2col関数が引っ掛かりました。
比較してみるとim2col関数にはSample数の次元がありますが、やってることはほとんど同じです。

結局、畳み込み演算も単なる行列積演算なんですね。

備考：

1. 畳み込みフィルターの大きさが(5,5)の場合

input = np.pad(input, [(2,2),(2,2),(0,0)], 'constant')  # (224,224,8)=>(228,228,8)
x = np.zeros((224,224,8,5,5))

for i, j in itertools.product(range(5), range(5)):
    x[:,:,:,i,j] = input[i:224+i,j:224+j]

2. stridesの大きさが(1,2)の場合

strides = (1,2)
input = np.pad(input, [(1,1),(1,1),(0,0)], 'constant')  # (224,224,8)=>(226,226,8)

for i, j in itertools.product(range(3), range(3)):
    x[:,:,:,i,j] = input[i:224+i:strides[0],j:224+j:strides[1]]

3. dilated convolutionの場合

input = np.pad(input, [(2,2),(2,2),(0,0)], 'constant')  # (224,224,8)=>(228,228,8)

for i, j in itertools.product(range(3), range(3)):
    x[:,:,:,i,j] = input[i*2:224+i*2,j*2:224+j*2]

目的

BigQueryにある緯度経度情報をTableauで可視化したい

問題

緯度経度情報（GEO系）をBigqueryでPOINTとかに変換してもTableauで文字列としてしか認識してくれない

解決方法

・python でshapeファイルに変換（geojsonとかでもなんでもいいと思う）

サンプルコード

import csv
import shapefile

w = shapefile.Writer("../data/point", shapeType=shapefile.POINT)

w.field("id", "C", "40")
w.field("date", "C", "40")
w.field("price", "C", "40")
w.field("flg", "C", "40")

with open("../data/___.csv") as text:
    reader = csv.DictReader(text)
    for row in reader:
        w.record(str(row['id']), str(row['date']), str(row['price']), str(row['flg']))
        w.point(float(row['longitude']), float(row['latitude']))

w.close()