np.random.seed(0)
X0 = sp.stats.norm(-2, 1).rvs(40)
X1 = sp.stats.norm(+2, 1).rvs(60)
X = np.hstack([X0, X1])[:, np.newaxis]
y0 = np.zeros(40)
y1 = np.ones(60)
y = np.hstack([y0, y1])

sns.distplot(X0, rug=True, kde=False, norm_hist=True, label="class 0")
sns.distplot(X1, rug=True, kde=False, norm_hist=True, label="class 1")
plt.legend()
plt.xlim(-6,6)
plt.show()

from sklearn.naive_bayes import GaussianNB
clf_norm = GaussianNB().fit(X, y)

clf_norm.classes_

clf_norm.class_count_

clf_norm.class_prior_

clf_norm.theta_, clf_norm.sigma_

xx = np.linspace(-6, 6, 100)
p0 = sp.stats.norm(clf_norm.theta_[0], clf_norm.sigma_[0]).pdf(xx)
p1 = sp.stats.norm(clf_norm.theta_[1], clf_norm.sigma_[1]).pdf(xx)                                                               
sns.distplot(X0, rug=True, kde=False, norm_hist=True, color="r", label="class 0 histogram")
sns.distplot(X1, rug=True, kde=False, norm_hist=True, color="b", label="class 1 histogram")
plt.plot(xx, p0, c="r", label="class 0 est. pdf")
plt.plot(xx, p1, c="b", label="class 1 est. pdf")
plt.legend()
plt.show()

x_new = -1
clf_norm.predict_proba([[x_new]])

px = sp.stats.norm(clf_norm.theta_, np.sqrt(clf_norm.sigma_)).pdf(x_new)
px

p = px.flatten() * clf_norm.class_prior_
p

clf_norm.class_prior_

p / p.sum()

"""
Explanation: 나이브 베이즈 분류 모형
나이브 베이즈 분류 모형(Naive Bayes classification model)은 대표적인 확률적 생성 모형이다. 
타겟 변수 $y$의 각 클래스 ${C_1,\cdots,C_K}$ 에 대한 독립 변수 $x$의 조건부 확률 분포 정보 $p(x \mid y = C_k)$ 를 사용하여 주어진 새로운 독립 변수 값 $x_{\text{new}}$에 대한 타켓 변수의 각 클래스의 조건부 확률 $p(y = C_k \mid x_{\text{new}})$ 를 추정한 후 가장 조건부 확률이 큰 클래스 $k$를 선택하는 방법이다.
조건부 확률의 계산
다음과 같이 베이즈 규칙을 사용하여 조건부 확률 $p(y = C_k \mid x_{\text{new}})$ 을 계산한다.
$$
P(y = C_k \mid x_{\text{new}}) = \dfrac{P(x_{\text{new}} \mid y = C_k)\; P(y = C_k)}{P(x_{\text{new}})}
$$
최종적으로는 각 클래스 $k$에 대한 확률을 비교하여 최고값을 계산하기만 하면 되므로 분모에 있는 주변 확률(marginal probability) ${P(x_{\text{new}})}$은 계산하지 않는다.
$$
P(y = C_k \mid x_{\text{new}}) \;\; \propto \;\; P(x_{\text{new}} \mid y = C_k) \; P(y = C_k)
$$
여기에서 사전 확률(prior) $P(y = C_k)$는 다음과 같이 쉽게 구할 수 있다.
$$
P(y = C_k) \approx \frac{\text{number of samples with }y = C_k}{\text{number of all samples}}
$$
$y$에 대한 $x$의 조건부 확률인 우도(likelihood)의 경우에는 일반적으로 정규 분포나 베르누이 분포와 같은 특정한 모형을 가정하여 다음과 같이 계산한다.

$P(x \mid y = C_k)$ 가 특정한 확률 분포 모형을 따른다고 가정한다.
트레이닝 데이터 ${x_1, \cdots, x_N}$을 사용하여 이 모형의 모수(parameter)를 구한다.
모수를 알고 있으므로 새로운 독립 변수 값 $x_{\text{new}}$이 어떤 값이 되더라도 $P(x_{\text{new}} \mid y = C_k)$ 를 계산할 수 있다.

우도 모형
우도의 모형으로 많이 사용하는 것은 다음과 같다.

베르누이 분포
$x$가 0 또는 1 값만을 가질 수 있다.
$x$가 1 이 될 확률은 고정되어 있다.
예: 동전을 던진 결과로 어느 동전을 던졌는지를 찾아내는 모형

$$ P(x_i \mid y = C_k) = \theta_k^x (1-\theta_k)^{(1-x_i)} $$

다항 분포 
$(x_1, \ldots, x_n)$ 이 0 또는 양의 정수
예: 주사위를 던진 결과로 어느 주사위를 던졌는지를 찾아내는 모형

$$ P(x_1, \ldots, x_n \mid y = C_k) = \prod_i \theta_k^{x_i}$$

가우시안 정규 분포
$x$가 실수로 특정한 값 근처
예: 시험 점수로 학생이 누구인지를 찾아내는 모형

$$ P(x_i \mid y = C_k) = \dfrac{1}{\sqrt{2\pi\sigma_k^2}} \exp \left(-\dfrac{(x_i-\mu_k)^2}{2\sigma_k^2}\right) $$
나이브 가정
독립 변수 $x$가 다차원(multi-dimensional) $x = (x_1, \ldots, x_n)$ 이면 위에서 사용한 우도 $P(x \mid y = C_k)$  는 원래 모든 $x_i$에 대한 결합 확률(joint probability) $P(x_1, \ldots, x_n \mid y = C_k)$ 을 사용해야 한다. 그러나 이러한 결합 확률은 실제로 입수하기 어렵기 때문에 모든 차원의 개별 독립 변수 요소들이 서로 독립(independent)이라는 가정을 흔히 사용한다. 이러한 가정을 나이브 가정(Naive assumption)이라고 한다.
나이브 가정하에서는 결합 확률이 개별 확률의 곱으로 나타난다.
$$
P(x_1, \ldots, x_n \mid y = C_k) = \prod_{i=1}^n P(x_i \mid y = C_k)
$$
$$
P(y = C_k \mid x_{\text{new}}) \;\; \propto \;\; \prod_{i=1}^n P(x_{\text{new},i} \mid y = C_k)\; P(y = C_k)
$$
Scikit-Learn에서 제공하는 나이브 베이즈 모형
Scikit-Learn의 naive_bayes 서브패키지에서는 다음과 같은 세가지 나이브 베이즈 모형 클래스를 제공한다.

BernoulliNB: 베르누이 분포 나이브 베이즈
MultinomialNB: 다항 분포 나이브 베이즈
GaussianNB: 가우시안 정규 분포 나이브 베이즈

이 클래스들은 다음과 같은 속성값 및 메서드를 가진다.


classes_: 공통


타겟 Y의 클래스(라벨)


class_count_: 공통


타겟 Y의 값이 특정한 클래스인 표본 데이터의 수


feature_count_: 베르누이 분포나 다항 분포 


타겟 Y의 값이 특정한 클래스이면서 독립 변수 X의 값이 1인 표본 데이터의 수 (베르누이 분포).


타겟 Y의 값이 특정한 클래스인 독립 변수 X의 값의 합 (다항 분포). 독립 변수 값이 1또는 0만 가지는 경우에는 표본 데이터의 수가 된다. 


class_prior_: 가우시안 정규 분포


타겟 Y의 무조건부 확률 분포 $ P(Y) $


class_log_prior_: 베르누이 분포나 다항 분포 


타겟 Y의 무조건부 확률 분포의 로그 $ \log P(Y) $


theta_, sigma_ : 가우시안 정규 분포


가우시안 정규 분포의 기댓값 $\mu$ 과 분산 $\sigma^2$


feature_log_prob_: 베르누이 분포나 다항 분포 


베르누이 분포 혹은 다항 분포의 모수 벡터의 로그
 $$ \log \theta = (\log \theta_1, \ldots, \log \theta_n) = \left( \log \dfrac{N_i}{N}, \ldots, \log \dfrac{N_n}{N} \right)$$


스무딩(smoothing)
 $$ \hat{\theta} = \frac{ N_{i} + \alpha}{N + \alpha n} $$ 


predict_proba(x_new) : 공통


조건부 확률 분포 $ P(Y \mid X_{\text{new}}) $


가우시안 정규 분포 나이브 베이즈 모형
End of explanation
"""


np.random.seed(0)
X = np.random.randint(2, size=(10, 4))
y = np.array([0,0,0,0,1,1,1,1,1,1])
print(X)
print(y)

from sklearn.naive_bayes import BernoulliNB
clf_bern = BernoulliNB().fit(X, y)

clf_bern.classes_

clf_bern.class_count_

np.exp(clf_bern.class_log_prior_)

fc = clf_bern.feature_count_
fc

fc / np.repeat(clf_bern.class_count_[:, np.newaxis], 4, axis=1)

theta = np.exp(clf_bern.feature_log_prob_)
theta

x_new = np.array([1, 1, 0, 0])

clf_bern.predict_proba([x_new])

p = ((theta**x_new)*(1-theta)**(1-x_new)).prod(axis=1)*np.exp(clf_bern.class_log_prior_)
p / p.sum()

x_new = np.array([0, 0, 1, 1])

clf_bern.predict_proba([x_new])

p = ((theta**x_new)*(1-theta)**(1-x_new)).prod(axis=1)*np.exp(clf_bern.class_log_prior_)
p / p.sum()

"""
Explanation: 베르누이 분포 나이브 베이즈 모형
베르누이 나이브 베이즈 모형에서는 타겟 변수뿐 아니라 독립 변수도 0 또는 1의 값을 가져야 한다. 예를 들어 전자우편과 같은 문서 내에 특정한 단어가 포함되어 있는지의 여부는 베르누이 확률 변수로 모형화할 수 있으므로 스팸 필터링에 사용할 수 있다.
End of explanation
"""


from sklearn.naive_bayes import MultinomialNB
clf_mult = MultinomialNB().fit(X, y)

clf_mult.classes_

clf_mult.class_count_

fc = clf_mult.feature_count_
fc

fc / np.repeat(fc.sum(axis=1)[:, np.newaxis], 4, axis=1)

clf_mult.alpha

(fc + clf_mult.alpha) / (np.repeat(fc.sum(axis=1)[:, np.newaxis], 4, axis=1) + clf_mult.alpha * X.shape[1])

theta = np.exp(clf_mult.feature_log_prob_)
theta

x_new = np.array([21, 35, 29, 14])
clf_mult.predict_proba([x_new])

p = (theta**x_new).prod(axis=1)*np.exp(clf_bern.class_log_prior_)
p / p.sum()

x_new = np.array([18, 24, 35, 24])
clf_mult.predict_proba([x_new])

"""
Explanation: 다항 분포 나이브 베이즈 모형
End of explanation
"""


from sklearn.datasets import fetch_20newsgroups
from sklearn.cross_validation import train_test_split
news = fetch_20newsgroups(subset="all")
X_train, X_test, y_train, y_test = train_test_split(news.data, news.target, test_size=0.1, random_state=1)

from sklearn.feature_extraction.text import TfidfVectorizer, HashingVectorizer, CountVectorizer
from sklearn.naive_bayes import MultinomialNB
from sklearn.pipeline import Pipeline

clf_1 = Pipeline([
            ('vect', CountVectorizer()), 
            ('clf', MultinomialNB()),
        ])
clf_2 = Pipeline([
            ('vect', TfidfVectorizer()),
            ('clf', MultinomialNB()),
        ])
clf_3 = Pipeline([
            ('vect', TfidfVectorizer(token_pattern=r"\b[a-z0-9_\-\.]+[a-z][a-z0-9_\-\.]+\b")),
            ('clf', MultinomialNB()),
        ])
clf_4 = Pipeline([
            ('vect', TfidfVectorizer(stop_words="english", 
                                     token_pattern=r"\b[a-z0-9_\-\.]+[a-z][a-z0-9_\-\.]+\b")),
            ('clf', MultinomialNB()),
        ])
clf_5 = Pipeline([
            ('vect', TfidfVectorizer(stop_words="english", 
                                     token_pattern=ur"\b[a-z0-9_\-\.]+[a-z][a-z0-9_\-\.]+\b")),
            ('clf', MultinomialNB(alpha=0.01)),
        ])

from sklearn.cross_validation import cross_val_score, KFold
from scipy.stats import sem

for i, clf in enumerate([clf_1, clf_2, clf_3, clf_4, clf_5]):
    scores = cross_val_score(clf, X_test, y_test, cv=5)
    print(("Model {0:d}: Mean score: {1:.3f} (+/-{2:.3f})").format(i, np.mean(scores), sem(scores)))

"""
Explanation: 예 1: 뉴스 그룹
End of explanation
"""


import codecs
def read_data(filename):
    with codecs.open(filename, encoding='utf-8', mode='r') as f:
        data = [line.split('\t') for line in f.read().splitlines()]
        data = data[1:]   # header 제외
    return data
train_data = read_data('/home/dockeruser/data/nsmc/ratings_train.txt')
test_data = read_data('/home/dockeruser/data/nsmc/ratings_test.txt')

X = zip(*train_data)[1]
y = zip(*train_data)[2]
y = np.array(y, dtype=int)

from sklearn.cross_validation import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X, y, train_size=10000, test_size=10000)
len(X_train), len(X_test)

from konlpy.utils import pprint
pprint((X[0], y[0]))

from konlpy.tag import Twitter
pos_tagger = Twitter()

%%time
def tokenize(doc):
    return ['/'.join(t) for t in pos_tagger.pos(doc, norm=True, stem=True)]

train_docs = [(tokenize(row[1]), row[2]) for row in train_data[:10000]]
tokens = [t for d in train_docs for t in d[0]]
import nltk
text = nltk.Text(tokens, name='NMSC')
mpl.rcParams["font.family"] = "NanumGothic"
plt.figure(figsize=(12,10))
text.plot(50)
plt.show()

"""
Explanation: 감성 분석 Sentiment Analysis

서울대 박은정님의 네이버 영화 감상평에 대한 감성 분석 예제
https://github.com/e9t/nsmc
https://www.lucypark.kr/slides/2015-pyconkr/

데이터 전처리
End of explanation
"""


from sklearn.feature_extraction.text import CountVectorizer
from sklearn.naive_bayes import MultinomialNB
from sklearn.pipeline import Pipeline
from sklearn.metrics import classification_report

clf_1 = Pipeline([
            ('vect', CountVectorizer()), 
            ('clf', MultinomialNB()),
        ])

%%time 
clf_1.fit(X_train, y_train)

pprint(list(clf_1.named_steps["vect"].vocabulary_)[:10])

%%time 
print(classification_report(y_test, clf_1.predict(X_test)))

"""
Explanation: CountVectorize 사용
End of explanation
"""


from sklearn.feature_extraction.text import TfidfVectorizer

clf_2 = Pipeline([
            ('vect', TfidfVectorizer()), 
            ('clf', MultinomialNB()),
        ])

%%time 
clf_2.fit(X_train, y_train)

%%time 
print(classification_report(y_test, clf_2.predict(X_test)))

"""
Explanation: TfidfVectorizer 사용
End of explanation
"""


from konlpy.tag import Twitter
pos_tagger = Twitter()

def tokenize_pos(doc):
    return ['/'.join(t) for t in pos_tagger.pos(doc, norm=True, stem=True)]

clf_3 = Pipeline([
            ('vect', CountVectorizer(tokenizer=tokenize_pos)), 
            ('clf', MultinomialNB()),
        ])

%%time 
clf_3.fit(X_train, y_train)

pprint(list(clf_3.named_steps["vect"].vocabulary_)[:10])

%%time 
print(classification_report(y_test, clf_3.predict(X_test), digits=4))

vect3 = clf_3.named_steps["vect"]
idx3 = np.array(np.argsort(vect3.transform(X_train).sum(axis=0)))[0]
voca3 = np.array(vect3.get_feature_names()).flatten()
pprint(voca3[idx3[-20:]].tolist())

"""
Explanation: 형태소 분석기 사용
End of explanation
"""


clf_4 = Pipeline([
            ('vect', TfidfVectorizer(tokenizer=tokenize_pos, ngram_range=(1,2))), 
            ('clf', MultinomialNB()),
        ])

%%time 
clf_4.fit(X_train, y_train)

%%time 
print(classification_report(y_test, clf_4.predict(X_test), digits=4))

"""
Explanation: 최적화
End of explanation
"""


