2018年12月29日土曜日

PythonでMCA(コレスポンデンス分析)

Pythonでコレポンをやる

macでやってる
$python3 plain_mca3.py cross_table_for_mca.csv


plain_mca3.pyの中身

import sys

# mca,pandasをインポート
import mca
import pandas as pd

# csvデータ読み込む
# index_col:行のインデックスに用いる列番号。 (デフォルト: None)
df = pd.read_csv(sys.argv[1],index_col=0)

# コレスポンデンス分析
# ncol = df.shape[1]
# Benzécri補正
# mca_ben = mca.MCA(df, ncols=ncol, benzecri=False, TOL=1e-8)
mca_ben = mca.MCA(df, benzecri=False, TOL=1e-8)


# Rowsのスコア(座標)を書き出す
result_row = pd.DataFrame(mca_ben.fs_r(N=2))
result_row.index = list(df.index)
print ("Rows:")
print(result_row)
print('\n', end='')

# Columnsのスコア(座標)を書き出す
result_col = pd.DataFrame(mca_ben.fs_c(N=2))
result_col.index = list(df.columns)
print ("Columns:")
print(result_col)
print('\n', end='')



# N(成分:固有値の数)の算出:表頭と表側の少ない方から1を引いた数にする
cnt_column = len(list(df.columns))
cnt_index = len(list(df.index))

if(cnt_column >= cnt_index) :
     cnt_eigenvalue = cnt_index-1
else :
    cnt_eigenvalue = cnt_column-1


# 固有値(eigenvalue)と寄与率(explained variance of eigen vectors)
data = {'value': pd.Series(mca_ben.L),
            'ratio': mca_ben.expl_var(greenacre=False, N=cnt_eigenvalue)}
columns = ['value', 'ratio']
table2 = pd.DataFrame(data=data, columns=columns).fillna(0)
table2.index += 1
table2.loc['Σ'] = table2.sum()
table2.index.name = 'Factor'
print ("Principal inertias(eigenvalues):")
print(table2)
print('\n', end='')



# 作図用ライブラリ
import matplotlib.pyplot as plt
import matplotlib
# import random as rnd #ラベル つけるときに使用

# Jupyterの中で表示したい場合は、プログラム初頭で、%matplotlib inlineとする。
# すると、インライン表示される(しなければ、別ウインドウが開く)。
# %matplotlib inline


# グラフのサイズを指定
plt.rcParams["figure.figsize"] = [7, 7]

fig, ax = plt.subplots()

# print(matplotlib.colors.cnames) #色の確認

# 表頭をプロット
result_col.plot(0, 1, kind='scatter', ax=ax, color='C0', s=20, marker="o")
for k, v in result_col.iterrows():
    ax.annotate(k, v)

# 表側をプロット
result_row.plot(0, 1, kind='scatter', ax=ax, color='#FFA500', s=20, marker='.')
for k, v in result_row.iterrows():
    ax.annotate(k, v)

# plt.rcParams['font.family'] = 'IPAexGothic' #全体のフォントを設定
# plt.rcParams['font.size'] = 12 #フォントサイズを設定 default : 12
# plt.rcParams['xtick.labelsize'] = 10 # 横軸のフォントサイズ
# plt.rcParams['ytick.labelsize'] = 10 # 縦軸のフォントサイズ
# matplotlib.font_manager._rebuild()

# X軸Y軸の目盛線とラベル
plt.axhline(0, color='gray')
plt.axvline(0, color='gray')
plt.xlabel('Factor 1')
plt.ylabel('Factor 2')

# 任意(図の設定)
# plt.figure(figsize=(4,4)) #図の設定
# plt.rcParams["font.size"] = 10 #なにかの指定


# 図を見てみる
plt.show()
----



cross_table_for_mca.csvの中身
----
,yamamoto,instant,twice,facecare,linestone
strong,80,30,50,90,70
beautiful,20,80,60,30,50
cute,20,90,70,10,60
clever,60,20,50,70,40
big,80,10,50,90,40
smart,60,30,20,70,40
charming,30,60,80,10,90
lovely,40,90,70,50,60
fresh,10,80,40,20,30
traditional,70,10,20,50,40
----



実行結果
----
Rows:
                    0         1
strong      -0.345675 -0.033748
beautiful    0.418075  0.078905
cute         0.584282 -0.040474
clever      -0.356028 -0.006797
big         -0.541367  0.010732
smart       -0.401739  0.149244
charming     0.390651 -0.354082
lovely       0.244705  0.090889
fresh        0.596563  0.282539
traditional -0.550051 -0.081841

Columns:

                  0         1
yamamoto  -0.514181 -0.033323
instant    0.624395  0.192199
twice      0.251629 -0.113170
facecare  -0.523083  0.164357
linestone  0.110477 -0.198570

Principal inertias(eigenvalues):
           value     ratio
Factor
1       0.197552  0.856996
2       0.023801  0.103250
3       0.007040  0.030540
4       0.002124  0.009214
Σ       0.230517  1.000000





クロス表の結果がこのcsvくらいのものであれば、上記コードできちんと結果が出る。
ところがモノによってはこれではきちんと出力されない場合がある。
その場合は、以下TOLの閾値を緩めておく必要がある。
mca_ben = mca.MCA(df, benzecri=False, TOL=1e-8)

で、最終的にいくつにしたかというと、この記述自体をなくしたか、もっと桁数増やしたと思う、たぶん。。。


2018年11月17日土曜日

exif情報の一括削除

ImageMagickでできる

フォルダ内のExif情報を削除したい場合はオプションの「-strip」で可能
該当フォルダに移動した後、一括操作の「mogrify」を使う
mogrify -strip *.jpg


イメージマジックでフォルダ内のExif情報を一括削除する
https://www.s-oj.com/ec-business/system/imagemagick%E3%82%92%E4%BD%BF%E3%81%A3%E3%81%A6%E3%82%B5%E3%83%96%E3%83%87%E3%82%A3%E3%83%AC%E3%82%AF%E3%83%88%E3%83%AA%E3%81%AEexif%E6%83%85%E5%A0%B1%E3%82%92%E4%B8%80%E6%8B%AC%E5%89%8A%E9%99%A4/

----

すべてのpng画像をjpgに変更
mogrify -format jpg *.png
(元画像は残る)

https://liginc.co.jp/394506

----

ImageMagick基本編:Mogrifyを用いた画像の一括処理について
http://shumilinux.blogspot.com/2016/06/imagemagick-mogrify.html

2018年10月13日土曜日

心理学の成り立ち


哲学の時代
自分の感覚(外部世界)と、知識は一致していた
ヒポクラテスの四体液説みたいな感じ
宗教の教えみたいな感じ(万物は神が創ったものだとか)

ところが、科学や技術が進歩してくると、どうもそうでもないらしいことが判明してきた

医学、物理学、天文学、科学いろいろ
コペルニクス、ガリレオ、ダーウィンとかとか

どうやら自分の感覚を疑わざるを得ないようだ

こうした知識・感覚と事実に矛盾があることから、
人の心の動きを明らかにしたいという学問的欲求、探究心が発生
心理学という学問が萌芽し始める

考えてるだけでは哲学の領域を出ない
科学的学問として成り立たせたい

哲学を科学化させたのはヴントの功績かも

数量化、客観化、機械化、解釈の排除、演繹
誰もが納得しやすい、わかりやすい、結果だけをみればOKというふうにしたい

客観性:再現可能性、反証可能性
信頼性:確実性(手堅さ)、監査可能性
内的妥当性:有意味性、真正性(迫真性ではない)
外的妥当性:転用可能性、一般化限定性


観察者であるか、行為者であるか
目の前で猫がエサを食べているのを
じっと見守って事実を伝えるか
エサを食べさせようとするか

前者がシャルコー
催眠で、体の症状だけでなく、心の症状が引き起こす事態があることを明確にした

後者がフロイト
催眠で治療しようと思った

(2009年3月メモ)



アメリカの成り立ち

中世抜きで発生した国

ピューリタン

  • カトリックとは反対で、ギリシャ・ローマ文明の知識
  • 建物もギリシャ・ローマ風、中世的なゴシックはなし

中世とは

  • 奴隷制度と騎士道
  • 中世ヨーロッパは1,000年かけて奴隷をなくした
  • 戦争のルールが変わった
  • 対等の敵と考えるか、悪者と考えるか
  • 機関銃という無差別殺人兵器に抵抗感を感じるかどうか

モンロー主義

  • 外の世界と関わらない
  • 自国産出の石油を無限に使えるから成り立つ、ただし60年代までだった
  • 70年代からは石油輸入国に(その後シェールガスが出てくるわけだが)
  • OPECや有色人種国にも気を使わないといけない

るつぼ

  • 日本は民族の吹き溜まり:うまく混ざってる
  • アメリカはその前段階なのかも:モザイク的に散らばってるだけ、混ざってはいない



自分の欲望に他律的

政府に賭博を禁止してほしい、あると自分がやっちゃうから

ピューリタンは他律的に抑える

  • オランダの家の、道に面した部屋の窓の広さ
  • 不道徳なことをしていないかを互いに監視
  • デスパレートな妻たちみたいな感じ

プロテスタントと違いカトリックは

  • 自分の罪を告解(懺悔)することで神の赦しが得られると考える
  • スペインの家は、中庭があって、内に開かれている

日本もカトリック的発想

  • 世間様と家の自治
  • 世間の縛りはあくまで家に対して、個人ではない
  • 家制度が活きてる時代は、個人が家に守られてたとも言える
  • だが今は家が守ってくれない、個人に直接攻め込まれる
  • (今の不寛容さは、カトリックからプロテスタントに移行してると見ることもできる)

(2009年1月メモ)


2018年7月1日日曜日

アサガオの数咲き作り

発芽したら
・本葉の6枚目と7枚目の間で摘芯
・双葉も摘む
・1,2,6枚目の脇芽も摘む

3,4,5枚目の脇芽を伸ばしていく

・6枚目の脇芽はグイグイ伸びてくるのでちょいちょい摘む
・子葉の両脇の芽も摘む


・3本ある子づるのそれぞれ2枚目と3枚目の間で摘芯

その後5~6日で孫づるが伸び出し、ここに蕾がつく、らしい

https://www.shuminoengei.jp/?m=pc&a=page_image_slideshow&target_report_id=7213&num=7
https://www.shuminoengei.jp/?m=pc&a=page_image_slideshow&target_report_id=7213&num=8


2018年6月24日日曜日

頻度流とベイズ流

ラプラスの悪魔
自然現象は「現在の位置と状態, そして現象の法則さえ分かれば未来も完璧に予測できる」と考えられていた。 だがしかし、未来を予測するのに必要な情報が足りないときに, 人間は何らかの信念に基づいて予測するしかない。

この主観的な信念を数値化したもの→主観確率

サイコロを振って5の目が出る確率→6分の1(どの目が出るかの条件は1〜6いずれも等しいだろうという信念が前提にある)

信念じゃないものに頼るなら、とりあえず振ってみる、何度も振ってみる。5の目が出た頻度を数え上げる。
過去の実績に基づいてるのが信念よりも確からしそうだし、誰の目から見ても変わらない実績なので信念よりもあてになりそう
→頻度流(頻度に基づいて確率を定義する)

 頻度流統計では最もな尤度のみを用いる→最尤法
 ベイズ流では事前分布+尤度

ただしこの考え方だと試行回数、すなわちサンプルサイズにいろいろ依存することになる。

サンプルサイズをあまり意識しないでよい。
事前確率はひとまず2分の1とでもおく。試行をちょっとずつ観察しながら5の目の出るもっともらしさ(尤度)の情報を得ながら事後確率を計算(更新)していく
→ベイズ流(主観確率に基づいて確率を定義する)


最尤法ではモデルが複雑な場合や、サンプルサイズよりも説明変数が多い場合などに無理があり、対応できない(サンプルサイズよりも説明変数が多いと最小二乗法、最尤法の計算が原理的に不可能になる)。

モデルが複雑になるとデータの最小単位がわかりにくい、サンプルサイズをカウントしづらい、解釈が難しい

厳密な実験計画が前提になっており、ベイズ更新、カルマンフィルタ(センサーなどから得られる情報をもとに機械の状態を把握し, リアルタイムで機械の動きを修正していくという工学の方法論で, 代表例としては人工衛星の軌道制御に応用される)のような後から追加情報が入ってくることを頻度流は想定してない

→頻度流の限界


反対に、ベイズ流の弱み

モデルが複雑で再現性が低い

事前分布をどう設定するかが初学者には難しい
(MCMCのようにたくさん事後確率を発生させるなら事前分布の有意味性が薄れる場合には問題視されにくい)



[教材] 今更だが, ベイズ統計とは何なのか.
http://ill-identified.hatenablog.com/entry/2017/03/17/025625

頻度論 vs. ベイズ統計(後半)
https://healthpolicyhealthecon.com/2015/12/28/bayesian-2/

コインで理解するベータ分布
http://r-tips.hatenablog.com/entry/beta-distribution

ベータ分布の確率密度、下側累積確率、上側累積確率のグラフを表示します。
https://keisan.casio.jp/exec/system/1161228837


600世帯での過去平均視聴率15%→見た人90人、見なかった人510人
(90,510)
この情報は、事前分布(1,1)+(89,509)の事後分布だったと考えることもできる

事前分布(1,1)にそれぞれ足して(91,511)で計算する

ベータ分布Beta(a, b)は「(a-1)回成功、(b-1)回失敗という情報が与えられたときの、二項分布のパラメータpの事後分布」と説明することができます(ただしa,bが非負整数のときのみ)。

ベイズ統計学とベータ分布 
http://www.creativ.xyz/beta-distribution-345