*Aya's Prelude*

大学院生のよしなごと

グレースケール画像の低ランク近似withSVDで遊んだ(追記:RGBもやってみた)

学校の授業で楽しそうな定理をやったので実際にやってみた.

定理

定理書こうと思ったけどはてなブログでLaTeX使うのめんどくさすぎる. やめた. 簡単にいうとある行列を特異値分解したときにそのUSVを使って低ランク k \leq nなる kで近似できるよって話. まじでLaTeX.

設定

実際に低ランクの行列でどこまででるのかやってみた. 用いる画像は「[甘美なる宝石]佐久間まゆ」の画像をグレースケールに変換したもの.サイズは800*640. これを以下のkランクまでの情報を持つ行列で近似していく.

結果

結果は以下のようになった.

  • 元の画像

f:id:remark_tzi:20191105234004j:plain
元画像

  • k = 600

f:id:remark_tzi:20191105233912j:plain
k = 600

  • k = 500

f:id:remark_tzi:20191105234039j:plain
k = 500

  • k = 400

f:id:remark_tzi:20191105234058j:plain
k = 400

  • k = 300

f:id:remark_tzi:20191105234112j:plain
k = 300

  • k = 200

f:id:remark_tzi:20191105234128j:plain
k = 200

  • k = 100

f:id:remark_tzi:20191105234145j:plain
k = 100

  • k = 50

f:id:remark_tzi:20191105234158j:plain
k = 50

気持ち

低ランクで結構絵になっててすげーという気持ちがある.

800640のデータが5050だけの情報であんなにちゃんと出るんだなぁという感動.

別に考察がしたかったりしたわけではないので楽しかったですって感想しか出ないな.

プログラム

環境はPython3. matlabでやればよかったけどまだPC変えてからいれてなかったのでまた今度.

# ライブラリのインポート
import numpy as np
from PIL import Image

# 画像の読み込みとグレースケール化
img = Image.open('images/mayu/mayu.jpg')
gray_img = img.convert('L')
gray_img.save('images/mayu/mayu_mono.jpg')

# 行列化
array = np.array(gray_img)

# 形状の確認
print(array.shape)

(800, 640)

# 落としたい次数(array.shapeで確認した数字以下)
k = 500

# 特異値分解
U,S,V = np.linalg.svd(array)
S = np.diag(S)

# 低ランク近似のための準備
U_k = U[:,:k]
S_k = S[:k,:k]
V_k = V[:k,:]

# 近似
array_k = np.dot(np.dot( U_k , S_k),V_k )

# 画像の保存
img_k = Image.fromarray(np.uint8(array_k))
img_k.save('images/mayu_mono_r%d.jpg' %k)

今後の展望

できれば今度はRGBとかでもやりたいなって思って色々考えたけどとりあえずRGBごとに特異値分解して元に戻すみたいなのしてもいいのかな？とか思った.

先生に色々相談したらRGBだとテンソルになるし高次元特異値分解ってのやってみなとのこと.

こっちもそういう定理あるのかな. 楽しそう.

追記

RGB版もできたら楽しいかなと思ってとりあえず行列ごとに分解して近似してみた.

設定

上の画像と同じもの. 但し, RGBのまま.

結果

kの設定も一緒！

  • 元の画像

f:id:remark_tzi:20191106002646j:plain
元の画像

  • k = 600

f:id:remark_tzi:20191106002702j:plain
k = 6

  • k = 500

f:id:remark_tzi:20191106002726j:plain
k = 500

  • k = 400

f:id:remark_tzi:20191106002739j:plain
k = 400

  • k = 300

f:id:remark_tzi:20191106002753j:plain
k = 300

  • k = 200

f:id:remark_tzi:20191106002805j:plain
k = 200

  • k = 100

f:id:remark_tzi:20191106002818j:plain
k = 100

  • k = 50

f:id:remark_tzi:20191106002833j:plain
k = 50

プログラム

# ライブラリのインポート
import numpy as np
from PIL import Image

# SVDと近似用の関数
def svd_and_approximation(array, rank):
    u,s,v = np.linalg.svd(array)
    s = np.diag(s)
    u_k = u[:,:rank]
    s_k = s[:rank,:rank]
    v_k = v[:rank,:]
    return np.asarray(np.dot(np.dot(u_k,s_k),v_k))

# 画像の読み込み
img = Image.open('images/mayu/mayu.jpg')

# 行列化
array = np.array(img)

# 形状の確認
print(array.shape)
print(array.size)
(800, 640, 3)
1536000

# 落としたい次数(array.shapeで確認した数字以下)
k = 600

# RGBごとに行列を分ける.
R = array[:,:,0]
G = array[:,:,1]
B = array[:,:,2]


array_R = svd_and_approximation(R,k)
array_G = svd_and_approximation(G,k)
array_B = svd_and_approximation(B,k)

rgb_k = np.asarray([array_R,array_G,array_B]).transpose(1,2,0)
 
# 画像の保存
img_k = Image.fromarray(np.uint8(rgb_k))
img_k.save('images/mayu_rgb/mayu_r%d.jpg' %k)

自然言語処理100本ノック Part1

mlp100本ノックやるぞ～～！！

やりたくね？という機運が高まったのでなんとなく始めた。
なんで始めたか書こうとしたらなんでかわからんくなった。そんな2019年春。
Part1ってなってるし続くと思うじゃん。
ブログにしようと思ったのがたまたまだから多分ブログには続かないんだな～～
サイトはこれ
ということでとりあえず10こやった。

Q0~Q10 準備運動

準備運動といってもきっとプログラムを普段書く人の準備運動だろうということで最初から気合い入れていくぞ～～！！

Q0 文字列の逆順

文字列"stressed"の文字を逆に（末尾から先頭に向かって）並べた文字列を得よ．

ぷよぐやむめっちゃやむ

a_0 = 'stressed'
print(a_0[::-1])

実行結果

'desserts'

いやまじでスライスとかいうの使う機会なさすぎて忘れてたから仕様を完全に忘れてめっちゃ時間かかった。ウケる。

Q1 「パタトクカシーー」

「パタトクカシーー」という文字列の1,3,5,7文字目を取り出して連結した文字列を得よ．


ぷよぐやむ

a_1 = 'パタトクカシーー'

print(a_1[::2])

上の実行結果

'パトカー'

ぷよぐやむ

print(a_1[1::2])

上の実行結果

'タクシー'

ここはさっきスライドの仕様思い出したからわかっちゃえば行けるなという気持ち

Q2 「パトカー」＋「タクシー」＝「パタトクカシーー」

「パトカー」＋「タクシー」の文字を先頭から交互に連結して文字列「パタトクカシーー」を得よ．

ぷよぐやむ

a_21 = 'パトカー'
a_22 = 'タクシー'
result = ''
for x,y in zip(a_21,a_22):
  result += x+y  
print(result)

実行結果

パタトクカシーー

うんちなので空っぽの文字列どうやってつくるのとか言い出す
空は青い

Q3 円周率

"Now I need a drink, alcoholic of course, after the heavy lectures involving quantum mechanics."という文を単語に分解し，各単語の（アルファベットの）文字数を先頭から出現順に並べたリストを作成せよ．
a_31 = "Now I need a drink, alcoholic of course, after the heavy lectures involving quantum mechanics."
a_32 = a_31.replace(',','').replace('.','').split( )
[len(x) for x in a_32]

実行結果

[3, 1, 4, 1, 5, 9, 2, 6, 5, 3, 5, 8, 9, 7, 9]

これはreplace()が出てこなくてなんだっっけ～～ってやった

Q4 元素記号

"Hi He Lied Because Boron Could Not Oxidize Fluorine. New Nations Might Also Sign Peace Security Clause. Arthur King Can."という文を単語に分解し，1, 5, 6, 7, 8, 9, 15, 16, 19番目の単語は先頭の1文字，それ以外の単語は先頭に2文字を取り出し，取り出した文字列から単語の位置（先頭から何番目の単語か）への連想配列（辞書型もしくはマップ型）を作成せよ．
a_41 = "Hi He Lied Because Boron Could Not Oxidize Fluorine. New Nations Might Also Sign Peace Security Clause. Arthur King Can."
a_42 = a_41.split( )

色々試してたやつ（だめだった）

one_frase = [0,4,5,6,7,8,14,18]
for i,(x,y) in enumerate(zip(a_42,one_frase)):
  print (i,x,y)

だめなやつの実行結果

0 Hi 0
1 He 4
2 Lied 5
3 Because 6
4 Boron 7
5 Could 8
6 Not 14
7 Oxidize 18

そっから集合にすればよくね？ってひらめいたやつ

one_frase = set(one_frase)
Element_symbol = {}
for i,x in enumerate(a_42):
  if i in one_frase:
    Element_symbol[i] = x[0]
  else :
    Element_symbol[i] = x[0:2]

よさげっぽくね？？プリントするぞい！！

print(Element_symbol)

実行結果

{0: 'H', 1: 'He', 2: 'Li', 3: 'Be', 4: 'B', 5: 'C', 6: 'N', 7: 'O', 8: 'F', 9: 'Ne', 10: 'Na', 11: 'Mi', 12: 'Al', 13: 'Si', 14: 'P', 15: 'Se', 16: 'Cl', 17: 'Ar', 18: 'K', 19: 'Ca'}

ここはなんか実験しながらこっちのがよさそ～～を繰り返してたらうまいこといってよかったなと思いました（小並感）
副産物としてenumerateとzipを組合す方法を覚えたので良いとする。

Q5 n-gram

与えられたシーケンス（文字列やリストなど）からn-gramを作る関数を作成せよ．この関数を用い，"I am an NLPer"という文から単語bi-gram，文字bi-gramを得よ．

文字列やリストなどなんて言われたもんだからビビっちゃってどうやって型確認するんだっけって試してたやつ

a = []
b = ''

type(a)

list

type(b)

str

試すの終わり
ちなみにクソほどわからんくてどうすればいいんだ～～～～～～～うわ～～～～～～～って言いながら書いた。絶対こんなに書かなくていいやつ。つーか自分で書いてていやだ～～～～～～～～！！！！って叫びながら書いてた。

ぷよぐやむ

def n_gram(n,word,sent_or_word):

  # 単語n-gramなのか文字n-gramなのかで前処理変えなきゃ
  if sent_or_word == True:
    word = word.split( )
  
  # 色々変数おいとこ
  start = 0
  results = []
  hokan = ''

  # n-gramします
  for j in range(len(word)):
    result = []
    for i in range(start,n+start):
        if j >= len(word)-n+1:
          break
        result.append(word[i])
    if j >= len(word)-n+1:
          break
    results.append(result)
    start += 1
  if sent_or_word == False:
    word_n_grams = []
    words = ''
    for item in results:
      for k in item:
        words += k
      word_n_grams.append(words)
      words = ''
    results = word_n_grams
  return results

n_gram(2,"I am the bone of my sword.",True)

[['I', 'am'],
 ['am', 'the'],
 ['the', 'bone'],
 ['bone', 'of'],
 ['of', 'my'],
 ['my', 'sword.']]

a_5 = 'I am an NLPer'

n_gram(2,a_5,False)

['I ', ' a', 'am', 'm ', ' a', 'an', 'n ', ' N', 'NL', 'LP', 'Pe', 'er']

先に感想書いちゃうくらいアホほど時間を使いました。
いやでも悩みながら書くのは楽しかったよ！

てか今気づいたけど3重ループじゃん草
なんとかしろよ

Q6 集合

"paraparaparadise"と"paragraph"に含まれる文字bi-gramの集合を，それぞれ, XとYとして求め，XとYの和集合，積集合，差集合を求めよ．さらに，'se'というbi-gramがXおよびYに含まれるかどうかを調べよ．

ぷよぐやむ

a_61 = "paraparaparadise"
a_62 = "paragraph"
X = set(n_gram(2,a_61,False))
Y = set(n_gram(2,a_62,False))

print(X)
print(Y)

実行結果

{'ar', 'ap', 'di', 'se', 'ad', 'ra', 'pa', 'is'}
{'ag', 'gr', 'ph', 'ar', 'ap', 'ra', 'pa'}

それぞれのもろもろの集合とったやつ

print(X | Y)
print(X & Y)
print(X - Y)
print(Y - X)

実行結果

{'ad', 'ag', 'ap', 'ar', 'di', 'gr', 'is', 'pa', 'ph', 'ra', 'se'}
{'ap', 'ar', 'pa', 'ra'}
{'ad', 'di', 'is', 'se'}
{'ag', 'gr', 'ph'}

ほんでこれが集合にあるかどうか確かめるやつ

if 'se' in X:
  print(True)
else :
  print(False)

True

if 'se' in Y:
  print(True)
else :
  print(False)

False

これ関数にしたら書くの一回で済んだのでは？？？？ウケてる

Q7 テンプレートによる文生成

引数x, y, zを受け取り「x時のyはz」という文字列を返す関数を実装せよ．さらに，x=12, y="気温", z=22.4として，実行結果を確認せよ．

関数名と変数名がカスで有名

def make_sentence(num_1,word,num_2):
  print(num_1,'時の',word,"は",num_2)

make_sentence(12,'気温',22.4)

実行結果

12 時の 気温 は 22.4

急に難易度が下がった気がしてしまったが実は俺の気づいてないギミックがあったのかもしれない。

Q8 暗号文

与えられた文字列の各文字を，以下の仕様で変換する関数cipherを実装せよ．

  ・ 英小文字ならば(219 - 文字コード)の文字に置換
  ・ その他の文字はそのまま出力

この関数を用い，英語のメッセージを暗号化・復号化せよ．

ここもう問題が何言ってるのかわからんかったわなんだよ文字コードって
めっちゃ調べたわ

def cipher(sentence):
  pre = ''
  for char in sentence:
    if char.isupper() == False:
      pre += chr(219 - ord(char))
    else :
      pre += char
  
  return pre

text = "I am the bone of my sword Steel is my body and fire is my blood I have created over a thousand blades"

Code = cipher(text)

Decryption = cipher(Code)

これはまじで忘れそうなので調べたサイトをきちんと乗っけるとする
小文字と大文字の区別のために使ったサイト
文字コードってなんだよ！！ってなったときのサイト

暗号化した結果

print(Code)

I zn gsv ylmv lu nb hdliw Sgvvo rh nb ylwb zmw uriv rh nb yollw I szev xivzgvw levi z gslfhzmw yozwvh

復号化した結果

print(Decryption)

I am the bone of my sword Steel is my body and fire is my blood I have created over a thousand blades

Q9 Typoglycemia

スペースで区切られた単語列に対して，各単語の先頭と末尾の文字は残し，それ以外の文字の順序をランダムに並び替えるプログラムを作成せよ．ただし，長さが４以下の単語は並び替えないこととする．適当な英語の文（例えば"I couldn't believe that I could actually understand what I was reading : the phenomenal power of the human mind ."）を与え，その実行結果を確認せよ．

最後の最後でまったくゆうことを聞かないタイプのポケモンなので例文は変えた

import random

def Typoglycemia(text):
  text = text.split( )

  results = ''
  
  for word in text:
    if len(word) <= 4:
      results += word + ' '
    else :
      sent = ''
      fi,la,mw = word[0],word[-1],list(word[1:-1])
      mws = random.sample(mw,len(mw))
      mwss = ''
      for k in mws:
        mwss += str(k)
      sent += fi + mwss + la + ' '
      results += sent
  return results

Typoglycemia(text)

'I am the bone of my sword Seetl is my body and fire is my bolod I have ceetard over a tunhosad belads '

書くの疲れたしもう言うことはない
しいていうならrandom.sampleの使い方かなぁ

今日のまとめ

思ったより考えたり書いたりするの楽しかった。
あと今日でめっちゃstr足せるじゃ～んってのとスライスだいぶ思い出せた。やってよかったな。

全体的な感想としては「へぇ～パズルかよ。」って感じだ。 楽しい。

GPMLの疑問点をほんとーに少し解消したつもりなので

21,41,48ページをなんとなく

  • ヤコビアンの理由
  • \Lambdaは対角行列？
  • なんで計画行列っていうの
  • 正規分布は情報量最小 の四点を解消したい。

まずヤコビアンだけど これ求めてるのが p(\mathbf{y})って気づいたので多分それで解消

積分するとこ変えたらそうなるなって感じ。

次に対角行列だけど やっぱり仮定に

"A^\mathrm{-1}を\lambdaの変換行列とする"

としないと厳しそう。

次に計画行列の由来

これはまじでググっただけでその人の感覚だったんだけどまぁあってるんじゃね？とかいう早朝5:30に見つけた記事をそのままはっつけるやつします（日本語やばすぎるけど明日起きて気が向いたら直すし記事もちゃんと読む。）（多分）

↓

計画行列（デザイン行列）とは何か：PRML編 - 北野坂備忘録

次に正規分布は情報量最小って話（エントロピー最大）

もうこれ読めば解決しそう（適当）↓

最大エントロピー法による正規分布の導出 - RigelのR言語メモであーる(R言語だけとは言っていない)

緑本の３章をpythonで勉強しました

データ解析のための統計モデリング入門３章

いわゆる緑本を勉強したのだが, 全部Rで書かれてたのでpythonのほうが好きかなーやっぱｗという気持ちで色々書いてみた.

モジュールのインポート

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import math
import matplotlib.style
import seaborn as sns

スタイルを変更したい気持ち

matplotlib.style.use("ggplot")

第3章 : 一般化線形モデル(GLM) -ポアソン回帰-

データの読み込み

data = pd.read_csv("data3a.csv")

データの確認

data.head()
y x f
0 6 8.31 C
1 6 9.44 C
2 6 9.50 C
3 12 9.07 C
4 10 10.16 C 
data["x"].head()

0     8.31
1     9.44
2     9.50
3     9.07
4    10.16
Name: x, dtype: float64
data["y"].head()

0     6
1     6
2     6
3    12
4    10
Name: y, dtype: int64
data["f"].head()

0    C
1    C
2    C
3    C
4    C
Name: f, dtype: object

3.2 観測されたデータの概要を調べる

data.describe()
y x
count 100.000000 100.000000
mean 7.830000 10.089100
std 2.624881 1.008049
min 2.000000 7.190000
25% 6.000000 9.427500
50% 8.000000 10.155000
75% 10.000000 10.685000
max 15.000000 12.400000 
fig,ax = plt.subplots() # fig は描画プロット, ax はサブプロット

X = data["x"]
y = data["y"]

ax.scatter(X,y, label = "x-y") # 散布図を出すときはscatter関数, label="hogehoge" とすることで凡例のラベルをつけられる
ax.set_xlabel("data_x") # x軸のラベル
ax.set_ylabel("data_y") # y軸のラベル
ax.set_title("data_x - data_y") # グラフの表題
ax.legend(loc = "best") # 凡例を表示する. loc="best"とすることで最も重なりのない位置に凡例を表示することができる.
plt.show()

f:id:remark_tzi:20190205040119p:plain

# 描画プロットとサブプロットの設定
fig,ax = plt.subplots() # fig は描画プロット, ax はサブプロット

## dataをf==Cとf==Tで分ける
data_C = data[data["f"]=="C"] # dataからf==Cのものだけをとってくる
data_T = data[data["f"]=="T"] # dataからf==Tのものだけをとってくる

## CとTごとにxとyに分ける
x_c = data_C["x"] 
y_c = data_C["y"]
x_t = data_T["x"]
y_t = data_T["y"]

## 散布図を出力する
ax.scatter(x_c,y_c,marker="v", label = "x_c-y_c") # marker = "v" とすると下三角
ax.scatter(x_t,y_t,marker="s", label = "x_t-y_t") # marker = "s" とすると正四角形, 他にも^や*やxなどがある
ax.set_xlabel("data_x") # x軸のラベル
ax.set_ylabel("data_y") # y軸のラベル
ax.set_title("data_x - data_y") # グラフの表題
ax.legend(loc = "best") # 凡例を表示する. loc="best"とすることで最も重なりのない位置に凡例を表示することができる.
plt.show()

f:id:remark_tzi:20190205040122p:plain

## f==cとf==Tの場合の箱ひげ図もみてみる
fig,ax = plt.subplots()

label_name = ["f = C","f = T"]
ax.boxplot((y_c,y_t), labels = label_name)

plt.show()

f:id:remark_tzi:20190205040126p:plain

## f==cとf==Tの場合の箱ひげ図もみてみる
fig,ax = plt.subplots()

#label_name = ["f = C","f = T"]
ax.boxplot(y_c,labels = "c")

plt.show()

f:id:remark_tzi:20190205040129p:plain

## f==cとf==Tの場合の箱ひげ図もみてみる
fig,ax = plt.subplots()

#label_name = ["f = C","f = T"]
ax.boxplot(y_c,labels = "ac")

plt.show()

---------------------------------------------------------------------------

ValueError                                Traceback (most recent call last)

<ipython-input-17-f0fe891b5ccf> in <module>
      3 
      4 #label_name = ["f = C","f = T"]
----> 5 ax.boxplot(y_c,labels = "ac")
      6 
      7 plt.show()
c:\users\nago\documents\python\studying\python3.6\lib\site-packages\matplotlib\__init__.py in inner(ax, data, *args, **kwargs)
   1808                         "the Matplotlib list!)" % (label_namer, func.__name__),
   1809                         RuntimeWarning, stacklevel=2)
-> 1810             return func(ax, *args, **kwargs)
   1811 
   1812         inner.__doc__ = _add_data_doc(inner.__doc__,
c:\users\nago\documents\python\studying\python3.6\lib\site-packages\matplotlib\axes\_axes.py in boxplot(self, x, notch, sym, vert, whis, positions, widths, patch_artist, bootstrap, usermedians, conf_intervals, meanline, showmeans, showcaps, showbox, showfliers, boxprops, labels, flierprops, medianprops, meanprops, capprops, whiskerprops, manage_xticks, autorange, zorder)
   3501 
   3502         bxpstats = cbook.boxplot_stats(x, whis=whis, bootstrap=bootstrap,
-> 3503                                        labels=labels, autorange=autorange)
   3504         if notch is None:
   3505             notch = rcParams['boxplot.notch']
c:\users\nago\documents\python\studying\python3.6\lib\site-packages\matplotlib\cbook\__init__.py in boxplot_stats(X, whis, bootstrap, labels, autorange)
   1179         labels = itertools.repeat(None)
   1180     elif len(labels) != ncols:
-> 1181         raise ValueError("Dimensions of labels and X must be compatible")
   1182 
   1183     input_whis = whis
ValueError: Dimensions of labels and X must be compatible

ここはなんでエラー吐くのかよくわからんからメモで残しといた. そのうち調べる（これは調べないですわ）

まあでも以上のことから

  • xが増加するに従ってyも増加しそう(よくわからない)
  • fはあまり影響しなさそう

とわかる

ポアソン回帰の統計モデル

ここで使うモジュールのインポート

## 統計的にはこのモジュールを使うと良いとされている
## というのも詳細がわかりやすい
import statsmodels.api as sm

GLMのフィッティング

x-yでポアソン回帰

# 植物の大きさ
x = sm.add_constant(data["x"])
# 飛ばした種子の数
y = data["y"]
# ポアソン回帰
model = sm.GLM(y, x, family=sm.families.Poisson()) # modelの作成
results = model.fit() 

results.summary()
Generalized Linear Model Regression Results
Dep. Variable: y No. Observations: 100
Model: GLM Df Residuals: 98
Model Family: Poisson Df Model: 1
Link Function: log Scale: 1.0000
Method: IRLS Log-Likelihood: -235.39
Date: Tue, 05 Feb 2019 Deviance: 84.993
Time: 03:27:42 Pearson chi2: 83.8
No. Iterations: 4 Covariance Type: nonrobust
coef std err z P>|z| [0.025 0.975]
const 1.2917 0.364 3.552 0.000 0.579 2.005
x 0.0757 0.036 2.125 0.034 0.006 0.145

f-yでポアソン回帰

## 施肥効果の変数
f = sm.add_constant(data["f"])
## 何もしてなければ0,なにかしていれば1に置き換える
f = f.replace({"C":0,"T":1})

## モデルの作成
model = sm.GLM(y, f, family=sm.families.Poisson())
results = model.fit()

## 詳細
results.summary()
Generalized Linear Model Regression Results
Dep. Variable: y No. Observations: 100
Model: GLM Df Residuals: 98
Model Family: Poisson Df Model: 1
Link Function: log Scale: 1.0000
Method: IRLS Log-Likelihood: -237.63
Date: Tue, 05 Feb 2019 Deviance: 89.475
Time: 03:27:44 Pearson chi2: 87.1
No. Iterations: 4 Covariance Type: nonrobust
coef std err z P>|z| [0.025 0.975]
const 2.0516 0.051 40.463 0.000 1.952 2.151
f 0.0128 0.071 0.179 0.858 -0.127 0.153

説明変数が２つ以上あるとき

## 説明変数をまとめる
Xs = data.drop("y",axis=1)
## objectoをダミー変数化
Xs = Xs.replace({"C":0,"T":1})
# 確率モデルをポアソン分布と仮定してモデルのフィッティング
Xs = sm.add_constant(Xs)
"""
# Xs = Xs.values # どうやらわざわざndarrayにしなくてもやってくれるらしい
"""
results = sm.GLM(endog=y, exog=Xs,family=sm.families.Poisson()).fit()
results.summary()
Generalized Linear Model Regression Results
Dep. Variable: y No. Observations: 100
Model: GLM Df Residuals: 97
Model Family: Poisson Df Model: 2
Link Function: log Scale: 1.0000
Method: IRLS Log-Likelihood: -235.29
Date: Tue, 05 Feb 2019 Deviance: 84.808
Time: 03:27:47 Pearson chi2: 83.8
No. Iterations: 4 Covariance Type: nonrobust
coef std err z P>|z| [0.025 0.975]
const 1.2631 0.370 3.417 0.001 0.539 1.988
x 0.0801 0.037 2.162 0.031 0.007 0.153
f -0.0320 0.074 -0.430 0.667 -0.178 0.114 
## 説明変数をまとめる
Xs = data.drop("y",axis=1)
## objectoをダミー変数化
Xs = Xs.replace({"C":0,"T":1})
# 確率モデルを正規分布と仮定してモデルのフィッティング
Xs = sm.add_constant(Xs)
"""
# Xs = Xs.values # どうやらわざわざndarrayにしなくてもやってくれるらしい
"""
results = sm.GLM(endog=y, exog=Xs,family=sm.families.Gaussian()).fit()
results.summary()
Generalized Linear Model Regression Results
Dep. Variable: y No. Observations: 100
Model: GLM Df Residuals: 97
Model Family: Gaussian Df Model: 2
Link Function: identity Scale: 6.6522
Method: IRLS Log-Likelihood: -235.12
Date: Tue, 05 Feb 2019 Deviance: 645.26
Time: 03:27:48 Pearson chi2: 645.
No. Iterations: 3 Covariance Type: nonrobust
coef std err z P>|z| [0.025 0.975]
const 1.6169 2.653 0.610 0.542 -3.582 6.816
x 0.6284 0.268 2.345 0.019 0.103 1.154
f -0.2538 0.537 -0.472 0.637 -1.307 0.800

正規分布とポアソン分布の線形回帰の違いを見てみる

poisson_const = 1.2631
poisson_f = -0.0320
poisson_consts = poisson_const + poisson_f
poisson_x =0.0801

j = np.arange(7.0, 12.1, 0.5)

poisson_link = [math.exp(poisson_consts + poisson_x * k) for k in j]

gauss_const = 1.6169
gauss_f = -0.2538
gauss_consts = gauss_const + gauss_f
gauss_x = 0.6284

gauss_link = [math.exp(gauss_consts + gauss_x * k) for k in j]
x_labels = [k for k in j]

gauss_id = [gauss_consts + gauss_x * k for k in j]

poisson_id = [poisson_consts + poisson_x * k for k in j]

fig, ax = plt.subplots()

#ax.scatter(x_labels,gauss_link, marker = "s", label = "Gauss_link")
ax.scatter(x_labels,poisson_link, marker = "*", label = "Poisson_link")
#ax.scatter(x_labels,gauss_id, marker = "v", label = "Gauss_id")
ax.scatter(x_labels,poisson_id, marker = "^", label = "Poisson_id")
ax.scatter(X,y,marker = "x", label = "data")
ax.set_xlabel("x_values")
ax.set_ylabel("lambda")
ax.legend(loc = "best") 
plt.show()

f:id:remark_tzi:20190205040114p:plain

タイトルに反して結局消してしまったガウス分布. 理由は簡単でこいつでやったのアホほど跳ね上がってほかがみえないからです. ごめんね！

スクレイピングを勉強する人 #色々な人

1. 実際にスクレイピングしてみた

1.1. はじめに

僕が無理言って友人にチュートリアル用に作ってもらったやつです。主にテキストしかとってこれないけどとりあえず慣れよう！

1.2. 目標

今回の目標はまゆのプロフィールをとってこよう！ってことですね！

f:id:remark_tzi:20180922064554p:plain
まゆのプロフィール

この画像の名前から下をとってくる。

1.3. ソースコードの説明

スクレイピングのソースコード書いてもらったので忘れないように説明を付けておく。

1.3.1. 今回使うモジュールのインポート

### スクレイピングするためのアレ

import requests as req
from lxml.html import *
import sys

ここは説明しなくてもいいよね。これから使うモジュールをインポートしただけ。どれをどこで使ってるかはようわからん。

1.3.2. Scrapyingするための関数の定義

classとかよくわからんけどとりあえずいろいろおまじないとかで今回はすます！

これはオブジェクト指向では必須ぽい。気持ちわかる。

### スクレイピングするためのアレ続き

class Scrapying:

    def __init__(self, url):
        self.res = req.get(url)
        self.result = {}

これは欲しい部分を探すやつです！ちなみに今欲しいのはこれ(青いところと灰色のところが対応してる)↓

f:id:remark_tzi:20180922064602p:plain
名前

f:id:remark_tzi:20180922064548p:plain
佐久間まゆ

### スクレイピングするためのアレ続き

    #tableタグを入手 #trまで行きたいってやつ
    def parse(self):
        if self.res.status_code != 200:
            print("Error: Page not found ")
            sys.exit()
        
        self.html = fromstring(self.res.text) #xpath使うときのおまじない
        self.html = self.html.xpath("//table//tr") #xpathの書き方は調べろ
  • //をつけると何もかもすっ飛ばして//の後に行く
  • //は間に何かないと使えない
  • 上の書き方するとリストになる。(xpathだと基本的にはリスト)

f:id:remark_tzi:20180922064545p:plain
みるとこ

今、ソースコードの紫のところがタグ。例で上げると

<div id="なんちゃら">

の<>で囲まれててかつ「紫＋何か」の構成のやつの紫部分。とりあえず名前と佐久間まゆの情報だけがほしいのでdivとかtableとかtrとかじゃなくてtdのところがほしいのでそこまで行こう！ってやつが

//table//tr

ってやつ。これは小技だけど一回tableとかtrとか経由することでうまいことリストが作れる。その辺は経験で覚えてこ。(まだよくわかってない)

### スクレイピングするためのアレ続き
               
    def after_parse(self):
        for html in self.html[1:]: #プロフィールってとこはいらないから1からにした
            title, content = html.xpath("td/text()") #td/text()ってのはtdタグに含まれるtextだけを取りに来る
  • 本来"タグA/タグB"は次のタグに行くだけど"タグA/text()"は別物
  • 今回は"名前"と"佐久間まゆ"みたいにtext部分が2つしかないので [名前,佐久間まゆ]っていうリストの中身をtitleとcontentにそれぞれぶち込んでる
### スクレイピングするためのアレ続き

            self.result[title] = content #resultは辞書配列なのでその辺は調べろ
            print("{}: {}".format(title, content)) #"{}: {}".format(A,B) のformatは{}にAとBをそれぞれぶち込むときの文法

ここは特にいうことない。コメントアウトよめばおけ。

### スクレイピングするためのアレ続き

if __name__ == '__main__': #おまじない
    s = Scrapying("https://www65.atwiki.jp/sakumamayu_wiki/")
    s.parse()
    s.after_parse()

ここも特になし。じゃあ切る必要あったかって言われるとよくわからん。

1.4. 出来上がったもの

f:id:remark_tzi:20180922064609p:plain
実行結果

って感じでプロフィールができた。嬉しい。作ったのは俺じゃないけどわかった！ってだけでうれしい。

1.5. Xpathの参考文献

これ見て頑張る。

1.6. 最後に

スクレイピングの練習でやってみただけなので特に悪いことをする予定はないです。(ちなみにスクレイピングで何かする予定も今はない)(そのうちなんかしたい)(しょうもないこととかしたいよね)モチベーション的にまゆが絡むとテンション上がるから選びました。モチベーションを作ってくれたこのサイトの運営者さんとまゆに感謝です。まゆすき。
みんなもそのサイトをみてまゆを知ろう→まゆの日記帳 - 佐久間まゆペディア @ wiki

スクレイピングしてみたいから勉強してる途中

1. スクレイピングの勉強をしてみた。

1.1. はじめに

はじめた理由は特にない。深夜テンション。楽しそうだった。

1.2. 今ある知識

プログラムとかスクレイピングとかの知識は0。しいていうならあんまりいっぱいアクセスしちゃいけないってことくらいしか知らない。なのでマジで0からスタートだから変なこと言ってたら教えてほしい。(優しいと嬉しい)

1.3. moduleなに使うか～～

とりあえずおすすめされたrequestsを使ってみようかなって思います。

1.4. その前に覚えとかないといけないこと！

スクレイピングするときはwebサーバーに負荷をかけないようにアクセスとアクセスの間には5秒位間隔をあけよう！！一番守らないといけないことらしい！！

1.5. やってみよー

とりあえずインストールしなきゃなのでインストールした。とりあえず参考にしたサイト(https://www.python.ambitious-engineer.com/archives/974)と同じようにやってみようと思ってやってみた。

1.5.1. 知らない言葉が出てきたからメモるやつ

  • ヘッダー
    ようわからんじゃ困るんだろうけどようわからん。なんならHTTPもようわからん。
  • get
    なんかとってくるのか？ようわからん。
  • post
    putとdeleteを補えるっぽい？ようわからん。
  • put
    編集するやつらしい。ようわからん。
  • エンコーディング
    ようわからんが文字コードがわかるっぽい。 これでとれるらしい。
requests.get(url).encoding
  • httpステータスコード
    これもようわからんけど友人にふわっときいたらサーバーがどうなってるか教えてくれるらしい。200なら正常っぽい？500とかだと鯖落ちしてるらしい。
  • レスポンスヘッダ
    ようわからんどころかなんもわからん。
  • byte形式
    いやなにこれ。

1.5.2. もはやわからんからちょっとわかるとこだけ説明したやつ

def get_encoding_from_headers(headers):
 
    content_type = headers.get('content-type')  # content-typeとは拡張子のこと(ほんとはちがう)(まじようわからん)
 
    if not content_type:                                      #なんもなかったらなんもないよってやつ？
        return None
 
    content_type, params = cgi.parse_header(content_type)      #ここはほんとなんもわからん
 

#こっからなんもわからん
    if 'charset' in params:
        return params['charset'].strip("'\"")
 
    if 'text' in content_type:       # 拡張子がtextだったら通す(あいまいな表現)
        return 'ISO-8859-1'

ここはなんもわからん。まあ気が向いたら頑張って調べるか。

1.6. そもそもわからんことが多すぎる件について

いやマジでわからんこと多すぎる。スクレイピング以前の問題説あるなってことで参考になりそうなサイトあったからそれを見ようと思う。
これ→https://review-of-my-life.blogspot.com/2017/10/python-scraping-beginners.html

はやくやりたいから頑張って読むぞー

1.7. 参考にしたサイトまとめ

はてなブログでもLaTeXのコマンドを使いたい!

1. はてなブログでも\LaTeXのコマンドを使いたい!

1.1. \LaTeXはいいぞ。

普段から何書くにしても\LaTeXを使っているのでどうしてもmdは使いにくいというか数式をかけなくてしんどい。ってことで調べてきた。

1.2. HTML編集できるやつだと・・・

下のやつをHTML編集の数式の一番上に入れると上手いこと行くらしい。やってないからしらんけど。まあそのうち試すかー（試さない）

<p style="display: none;">[tex: ]</p>

1.3. 肝心な方

markdownの方はどうなのかっていうと

[tex:ここに数式]

ってやるとうまいこと行く。これは試した。だけどこれは\LaTeXでいうとこの$数式$と同じ扱いだから

[tex:\sum_{k=1}^{n}k]

とかするとおっきくならない。実際これはこうなる \sum_{k=1}^{n}k
いちいち

[tex:\displaystyle\sum_{k=1}^{n}k]

ってしないとだめ。だるすぎか？なんかこうガーーって

\begin{align}
\sum_{k=1}^{n}k
\end{align}

いつも通り書く方法はないだろうか・・・もうちょっと調べてみようと思う。

1.4. \LaTeXを使う上でとっても強いサイト

いっつも参考にしているサイトというかチートシートというかを置いときます。いやまじでこれ便利。