단회귀분석

12652 단어 Python3
"Stan 및 R의 Base 통계 모델링"Chapter 07

비선형 관계


2차 커브 예


한 가정의 여름부터 겨울까지 매일 에어컨 소모 전력 YkWh와 실외의 평균 기온 X(C)의 기록
d = pd.read_csv("input/data-aircon.txt")

sns.set_style("whitegrid")
sns.jointplot(x="X",y="Y",data=d)

모형식


$ Y[n]\sim Normal(a+b(X[n] - x_0)^2,\sigma_Y)$　　　$ n=1,...,N $
에어컨 소모 전력 감소 $20\,^{circ}\mathrm{C} 달러, 이 가정의 편안한 온도 $x0달러(X-x 0)로 가정^2달러 증가
stanmodel = StanModel(file="model/model7-3.stan")

Area_new = np.linspace(-3,32,60)
data_ = {'N':len(d),'X':d['X'],'Y':d['Y'],'N_new':60,'X_new':Area_new}
fit1 = stanmodel.sampling(data=data_,n_jobs=-1,seed=1234)

col =  np.linspace(-3,32,60)
print (len(col))
df2 = pd.DataFrame(fit1['y_new'])
df2.columns = col

qua = [0.1, 0.25, 0.50, 0.75, 0.9]
d_est = pd.DataFrame()

for i in np.arange(len(df2.columns)):
    for qu in qua:
        d_est[qu] = df2.quantile(qu)

x = d_est.index
y1 = d_est[0.1].values
y2 = d_est[0.25].values
y3 = d_est[0.5].values
y4 = d_est[0.75].values
y5 = d_est[0.9].values


plt.fill_between(x,y1,y5,facecolor='blue',alpha=0.2)
plt.fill_between(x,y2,y4,facecolor='blue',alpha=0.5)
plt.plot(x,y3,'k-')
plt.scatter(d["X"],d["Y"],c='k',alpha=0.3)
plt.xlabel("X")
plt.ylabel("Y")
sns.set_style('whitegrid')
plt.title('Fig7_5')

예제 시간 시퀀스 데이터


환자에게 약물을 투여할 때 약을 복용한 후 시간타임(hour)과 약의 혈액농도 Y(mg/mL)를 거친 데이터
d = pd.read_csv("input/data-conc.txt")

x = d['Time'].values
y = d['Y'].values
plt.plot(x,y,marker='o')
plt.xlim(0,25)
plt.ylim(0,15)
plt.xlabel('Time(hour)')
plt.ylabel('Y')

모형식


$ Y[n]\sim Normal(1-exp(-b Time[t])),\sigma_Y) $　　　　　 $ t=1,,,,T $
배경 지식에 근거하여 곡선 파라미터에 제한을 가했다
$ real a;$
$real b; $
stanmodel = StanModel(file="model/model7-4.stan")

Time_new = np.linspace(0,24,60)
data_ = {'T':len(d),'Time':d['Time'],'Y':d['Y'],'T_new':60,'Time_new':Time_new}
fit1 = stanmodel.sampling(data=data_,n_jobs=-1,seed=1234)
ms1 = fit1.extract()

col =  np.linspace(0,24,60)
print (len(col))
df2 = pd.DataFrame(fit1['y_new'])
df2.columns = col

qua = [0.1, 0.25, 0.50, 0.75, 0.9]
d_est = pd.DataFrame()

for i in np.arange(len(df2.columns)):
    for qu in qua:
        d_est[qu] = df2.quantile(qu)

x = d_est.index
y1 = d_est[0.1].values
y2 = d_est[0.25].values
y3 = d_est[0.5].values
y4 = d_est[0.75].values
y5 = d_est[0.9].values


plt.fill_between(x,y1,y5,facecolor='blue',alpha=0.2)
plt.fill_between(x,y2,y4,facecolor='blue',alpha=0.5)
plt.plot(x,y3,'k-')
plt.scatter(d["Time"],d["Y"],c='k',alpha=0.3)
plt.xlabel("Time(hour)")
plt.ylabel("Y")
sns.set_style('whitegrid')
plt.title('Fig7_5')

시간 시퀀스 데이터의 비선형 함수 응용 예시


모델 1: 지수 함수 시계추의 진동 폭은 공기 저항 등에 따라 지수 함수로 감소한다
Model2: S형 증가 방식의 곡선
Model 3: 초기 증가, 커브 감소
y1 = lambda x: 2*np.exp(-x)
y2 = lambda x: 1.8/(1+50*np.exp(-2*x))
y3 = lambda x: 8*(np.exp(-x) - np.exp(-2*x))

x = np.linspace(0,5,100)
plt.plot(x,y1(x),'k-',label='1')
plt.plot(x,y2(x),'k--',label='2')
plt.plot(x,y3(x),'k.',label='3')
plt.legend()
plt.xlabel("Time")
plt.ylabel("Y")
plt.title('Fig7_7')

다중 공통 선형


multiplicollinearity: 재회귀분석에서 변수 간의 관련성이 매우 높다는 것을 설명하기 때문에 회귀계수는 수렴하지 않는다
변수의 배경 지식이 명확하게 관련되어 있는 상황에서 원칙적으로 한 쪽을 버리다

연락하다


confounding 모델 외에 응답 변수와 설명 변수에 영향을 주는 변수가 존재합니다
허구의 초등학생 50m 달리기 데이터, 평균 초속 $Y(m/sec)$, 변수는 체중 Weight(kg)
d = pd.read_csv("input/data-50m.txt")
x = d['Weight'].values
y = d['Y'].values
c = d['Age']

plt.scatter(x,y,c=c,edgecolors='k')
plt.xlabel("Weight")
plt.ylabel("Y")
sns.set_style('whitegrid')
plt.title('Fig7_8')


sns.pairplot(d, hue='Age',
             vars=['Weight', 'Y'])

모형식


$\mu_{weight}[n] = c_1 + c_2 Age[n] $
$ Weight[n]\sim Normal(\mu_{weight}[n],\sigma_W) $
$\mu_Y [n] = b_1 + b_2 Age[n] + b_3 Weight[n] $
$ Y[n]\sim Normal (\mu_Y [n],\sigma_Y ) $
변수 간의 인과 관계를 탐색하는 해석을 경로 해석이라고 한다
stanmodel = StanModel(file="model/model7-5.stan")
data_ = {'N':len(d),'Weight':d['Weight'],'Y':d['Y'],'Age':d['Age']}
fit1 = stanmodel.sampling(data=data_,n_jobs=-1,seed=1234)
ms1 = fit1.extract()

설명 변수가 너무 많음


설명 변수에 노이즈 포함

d = pd.read_csv("input/data-salary.txt")

모형식


$ X[n]\sim Normal(x_{true}[n],2.5 $
$ Y[n]\sim Normal(a + b X_{true}[n],\sigma_Y) $
stanmodel = StanModel(file="model/model7-6.stan")

data = {"N":len(d),"X":d["X"],"Y":d["Y"]}
fit = stanmodel.sampling(data=data,n_jobs=-1,seed=1234)

col = d['X'].values
print (len(col))
df2 = pd.DataFrame(fit['x_true'])
df2.columns = col

qua = [0.1, 0.25, 0.50, 0.75, 0.9]
d_est = pd.DataFrame()

for i in np.arange(len(df2.columns)):
    for qu in qua:
        d_est[qu] = df2.quantile(qu)

x = d_est.index
y1 = d_est[0.1].values
y2 = d_est[0.25].values
y3 = d_est[0.5].values
y4 = d_est[0.75].values
y5 = d_est[0.9].values

fig = plt.figure(figsize=matplotlib.figure.figaspect(1))

plt.fill_between(x,y1,y5,facecolor='blue',alpha=0.1)
plt.fill_between(x,y2,y4,facecolor='blue',alpha=0.5)
plt.plot(x,y3,'k-')
plt.scatter(x,x,c='b')
plt.xlabel("X")
plt.ylabel("X")
sns.set_style('whitegrid')
plt.title('Fig7_6')

마감하다


모형식


마감이 없는 경우.
$ Y[n]\sim Normal(\mu,\sigma_Y) $ 　　　$ n = 1,.., N_obs $
마감된 경우가 있어요.
$y[n]\sim 정규(\mu,\sigma Y)$하지만 $y[n]절단된 경우 측정값이 $y$ Prob[y$ =\int_{-\infty}^L\frac{1}{\sqrt{2\pi\sigma}} exp[-\frac{1}{2}(\frac{y -\mu}{\sigma})^2] d_y $
$ =\int_{-\infty}^{\frac{L-\mu}{\sigma}}\frac{1}{\sqrt{2\pi}} exp[-\frac{1}{2}z^2] d_z =\int_{-\infty}^{\frac{L-\mu}{\sigma}}\phi(z)d_z $
$ =\Phi (\frac{L -\mu}{\sigma}) $
$\phi 달러는 표준 정적 분포이고 $\Phi 달러는 누적 분포 함수cumulative distribution function(CDF)
d = pd.read_csv("input/data-protein.txt")

stanmodel = StanModel(file="model/model7-7.stan")

Y_obs = d[~d['Y'].str.contains('<')].astype('float')
L_ = d[d['Y'].str.contains('<')]
L = L_['Y'].str.replace('<','').astype('float')[0]
data = {"N_obs":len(Y_obs),"N_cens":len(L_),"Y_obs":Y_obs['Y'],"L":L}
fit = stanmodel.sampling(data=data,n_jobs=-1,seed=123)

편차 값:outlier


보기 드문 큰 값을 생성하는 메커니즘을 가정하고 편차 값을 포함하여 분석하는 방법

모형식


밑단 길이의 밑줄 분포와 Student의 t 분포를 사용합니다
$ Y[n]\sim Cauchy(a+bX[n],\sigma) $ 　　　$ n=1,2,...,N $
혼합 정적 분포와 Zero-Inflated Poisson 분포 등을 사용합니다.
d = pd.read_csv("input/data-outlier.txt")

plt.scatter(x=d['X'],y=d['Y'])
plt.xlabel("X")
plt.ylabel("Y")
sns.set_style('whitegrid')
plt.title('Fig7_11')


import hashlib

def StanModel_cache(model_code, model_name=None, **kwargs):
    """Use just as you would `stan`"""
    code_hash = hashlib.md5(model_code.encode('utf-8')).hexdigest()
    if model_name is None:
        cache_fn = path + 'stan_pkl/cached-model-{}.pkl'.format(code_hash)
    else:
        cache_fn = path + 'stan_pkl/cached-{}-{}.pkl'.format(model_name, code_hash)

    try:
        sm = pickle.load(open(cache_fn, 'rb'))
    except: 
        sm = StanModel(file=model_code)
        with open(cache_fn, 'wb') as f:
            pickle.dump(sm, f)
    else: 
        print("Using cached StanModel")
    return sm

model = StanModel_cache(model_code='model/model4-4.stan')


model_code='model/model4-4.stan'
code_hash = hashlib.md5(model_code.encode('utf-8')).hexdigest()

with open(path + 'stan_pkl/cached-model-{}.pkl'.format(code_hash), 'rb') as f:
    model = pickle.load(f) 

data = pd.read_csv("input/data-outlier.txt")

stan_data = data.to_dict('list')

X_new = np.arange(0, 11, 1)
stan_data.update({'N':len(data), 'X_new':X_new, 'N_new':len(X_new)})

fit_nuts = model.sampling(data = stan_data, seed=1234, warmup=1000, chains=4)

col = np.arange(11)
print (len(col))
df2 = pd.DataFrame(fit_nuts['y_new'])
df2.columns = col

qua = [0.025, 0.25, 0.50, 0.75, 0.975]
d_est = pd.DataFrame()

for i in np.arange(N_X):
    for qu in qua:
        d_est[qu] = df2.quantile(qu)

x = d_est.index
y1 = d_est[0.025].values
y2 = d_est[0.25].values
y3 = d_est[0.5].values
y4 = d_est[0.75].values
y5 = d_est[0.975].values


plt.fill_between(x,y1,y5,facecolor='blue',alpha=0.1)
plt.fill_between(x,y2,y4,facecolor='blue',alpha=0.5)
plt.plot(x,y3,'k-')
plt.scatter(d["X"],d["Y"],c='b')

plt.xlabel("X")
plt.ylabel("Y")
sns.set_style('whitegrid')
plt.title('Fig7_12')
경사가 편리치에 의해 견인되다

stanmodel = StanModel(file="model/model7-9.stan")

data = {"N":len(d),"X":d["X"],"Y":d["Y"],"N_new":len(np.arange(11)),"X_new":np.arange(11)}
fit = stanmodel.sampling(data=data,n_jobs=-1,seed=1234)

qua = [0.025, 0.25, 0.50, 0.75, 0.975]
d_est = pd.DataFrame()

for i in np.arange(N_X):
    for qu in qua:
        d_est[qu] = df2.quantile(qu)

x = d_est.index
y1 = d_est[0.025].values
y2 = d_est[0.25].values
y3 = d_est[0.5].values
y4 = d_est[0.75].values
y5 = d_est[0.975].values


plt.fill_between(x,y1,y5,facecolor='blue',alpha=0.1)
plt.fill_between(x,y2,y4,facecolor='blue',alpha=0.5)
plt.plot(x,y3,'k-')
plt.scatter(d["X"],d["Y"],c='b')

plt.xlabel("X")
plt.ylabel("Y")
sns.set_style('whitegrid')
plt.title('Fig7_12')
소음 부분은 해양 분포의 모형을 사용한다
경사도가 편차값에 끌리지 않았다
그러나 어느 것이 옳은 것은 아니다
Next.js 인증을 사용하여 Next.js 웹사이트에 간단한 인증 구현
정적 분석 도구 + FireLinker

좋은 웹페이지 즐겨찾기

개발자 우수 사이트 수집

개발자가 알아야 할 필수 사이트 100선 추천 우리는 당신을 위해 100개의 자주 사용하는 개발자 학습 사이트를 정리했습니다
best100-homepage