역문제 vs 순문제


역문제: 역문제란 $y=ax$와 같이 $y$을 알고, 원인 $x$을 추정하는 형태의 문제를 역문제라고함.

순문제: 순문제는 역문제와 반대로 x에서 y을 예측하는 문제를 순문제라고 함.

 

정칙성과 행렬


정칙행렬(Nonsingular matrix, regular matrix): 정칙행렬은 역행렬이 존재하는 정방행렬을 말한다. 

특이행렬(Singular matrix): 정칙행렬과 반대로 역행렬이 존재하지 않는 경우를 말한다.

 

 

연립일차방정식의 해법


$y=ax$와 같이 연립방정식의 해가 있는 경우는 A을 구할 수 있다. 즉, A가 정칙행렬인 경우에는 x을 구할 수 있다는 말이다.

연립일차방정식의 풀이는 정칙행렬인 경우 1) 변수소거법, 2)가우스 요르단 소거법, 3)역행렬을 계산으로 풀 수 있다.

 

1) 변수 소거법의 핵심은 $(I|S)(\frac{x}{-1})=o$ 을 만드는 것이다. 예를 들어, $A=\begin{bmatrix}
2 &  3& 3\\ 
3 &  4& 2\\ 
-2 &  -2& 3
\end{bmatrix}$, $y=\begin{bmatrix}
9\\ 
9\\ 
2
\end{bmatrix}$ 인 경우

- Step 1) $=O$꼴로 만든다.  \begin{bmatrix}
2 &  3&  3& 9\\ 
3 &  4&  2& 9\\ 
-2 &  -2&  3& 2 
\end{bmatrix} \begin{bmatrix}
x_{1}\\ 
x_{2}\\ 
x_{3}\\ 
-1
\end{bmatrix}=\begin{bmatrix}
0\\ 
0\\ 
0
\end{bmatrix}

- Step 2) A을 단위 행렬 로 만든다. 아래와 같이 대각행렬이 I인 행렬을 만든다.

=> 요약하면 $( I | S) (\frac{x}{-1})=o$의 꼴로만드는 것이 골자이다. (S은 상수들, |은 행렬의 이어붙이기)

다음은 numpy을 이용한 변수소거법의 예시이다.

import numpy as np

A = np.array([[2,3,3],[3,4,2],[-2,-2,3]], dtype=np.float32) # (3,3)
y = np.array([[9,9,2]], dtype=np.float32).T  # 열벡터 (3,1)
x = np.array([0,0,0], dtype=np.float32)

Ay = np.concatenate([A, y], axis=1)
print(Ay)
[[ 2.  3.  3.  9.]
 [ 3.  4.  2.  9.]
 [-2. -2.  3.  2.]]
 
 
Ay[0, :] = 1/2 * Ay[0, :]
print(Ay)
[[ 1.   1.5  1.5  4.5]
 [ 3.   4.   2.   9. ]
 [-2.  -2.   3.   2. ]]
 

Ay[1, :] = Ay[1, :] + -3 *Ay[0, :]
print(Ay)
[[ 1.   1.5  1.5  4.5]
 [ 0.  -0.5 -2.5 -4.5]
 [-2.  -2.   3.   2. ]]
 
Ay[1, :] = -2*Ay[1, :]
print(Ay)   # 행에 (-2배)
[[ 1.   1.5  1.5  4.5]
 [-0.   1.   5.   9. ]
 [-2.  -2.   3.   2. ]]
 
# 3행의 첫원소도 0으로 바꿔줌
Ay[2, :] = Ay[2, :] +2*Ay[0, :]
print(Ay)
[[ 1.   1.5  1.5  4.5]
 [-0.   1.   5.   9. ]
 [ 0.   1.   6.  11. ]]
 
# 3행의 두번째 원소도 0으로 바꿔줌
Ay[2, :] = Ay[2, :] - Ay[1, :]
print(Ay)
[[ 1.   1.5  1.5  4.5]
 [-0.   1.   5.   9. ]
 [ 0.   0.   1.   2. ]]
 
# 2행의 마지막원소 [5]을 [0]으로
Ay[1, :] = Ay[1, :] - 5*Ay[2, :]
print(Ay)
[[ 1.   1.5  1.5  4.5]
 [-0.   1.   0.  -1. ]
 [ 0.   0.   1.   2. ]]
 
# 1행의 1원소 [1] 제외하고 모두 [0]으로
Ay[0, :] = Ay[0,:] - 1.5 * Ay[1, :] - 1.5*Ay[2,:]
print(Ay)
[[ 1.  0.  0.  3.]
 [-0.  1.  0. -1.]
 [ 0.  0.  1.  2.]]
 
x = Ay[:, -1]
print(x)

[ 3. -1.  2.]
print(A@x)
[9. 9. 2.]

 

피보팅(Pivoting):  A의 꼴을 I로 만드는 가우스-요르단 소거법으로 연립방정식을 풀 때, 각 행의 앞의 원소가 0이어서 1만들수 없는 경우. $\begin{bmatrix}
0 & 1 & 6 & 11\\ 
3 & 2 & 5 & 7\\ 
5 & 2 & 1 & 2
\end{bmatrix} \rightarrow \begin{bmatrix}
3 & 2 & 5 & 7\\ 
0 & 1 & 6 & 11\\ 
5 & 2 & 1 & 2
\end{bmatrix}$

핵(Kernel): 주어진 A에 의해 $Ax=o$으로 이동해오는것과 같은 x의 집합을 A의 kernel이라고함. 표기: $Ker A$

상(Image): 주어진 A에 대해 x을 여러모로 움직인 경우에 A로 옮기는 y=Ax의 집합을 A의 상. 원래공간은 A을 이용해 옮긴 영역을 말함. 표기 $Im A$

 

특이행렬.


연립방정식에서 주어진 단서가 일치하는 경우. 예를 들어, 찾고자하는 x는 4개의 원소이고, 알고있는 y, A의 원소도 4개(4x4)이다. 그런데, 사실A의 단서인 열들이 겉보기에는 단서가 4개이지만, 사실 1개는 중복이어서 실질적으로 3개인 것을 의미한다.

이러한 단서의 실질적인 갯수를 랭크(rank) 라고한다.

A = np.array([[2,3,3],[3,4,2],[-2,-2,3]], dtype=np.float32) # (3,3)
rank_A = np.linalg.matrix_rank(A)
print(rank_A)

# 3

일대일 사상: "같은 결과 y결과가 나오는 원인 x가 유일한가 ?"

위로의 사상: "어떤 결과 y가 나와도, 그 원인 x가 무조건 존재하는가"

반응형

백터표기

1. 벡터는 기본적으로 종벡터(세로로 길게쓴 벡터)를 사용한다.

이유: 벡터는 주로 행렬과의 연산과 함께 쓰는데, 어순을 지키기 위함이다. 예를 들어, 행렬의 곱 $A\textbf{x}$을 어순을 바꿔쓰면, $\textbf{x}A$으로 부자연스럽기 때문이다.

2. 벡터는 주로 $\textbf{x, v, e}$와 같이 두꺼운 글씨로 표기한다.

이유: 벡터와 숫자를 확실히 구분하려고 한다. 기본이지만 이를 익혀두는 것이 좋다.

3. 숫자의 나열로서의 벡터는 $x$, 공간의미를 강조하고자할 때는 $\vec{x}$로 표기한다.

4. 좌표는 '기저'에 기인하지만, 극좌표의 경우도 있을 수 있다.

5. 아핀공간이란? 선형 공간과 달리, 공간에서의 원점이 없는 공간

6. 기저의 조건

 - 어떤 벡터 $\vec{v}$라도  $\vec{v} = x_{1}\vec{e_{1}}+...+x_{n}\vec{e_{n}}$로 표현할 수 있다.

 - 위와 같이 표현하는 방법으로 딱 한가지로만 표현할 수 있다.

7. 차원

 = 기저벡터의 개수 = 좌표의 성분의 수 

 

 

행렬: 수를 직사각 형태로 나타낸 것

1. 행렬의 표기: 알파벳 대문자로 표기하며, 두꺼운 글씨가아닌 단순한 대문자로 쓰는 것이 일반적.

2. 행렬의 곱

 - 행렬과 벡터의 곱은 벡터

 - 행렬의 열의 수(컬럼 수)가 차원수, 행 수(로우 수, 높이)가 출력의 차원수

 - 행렬의 A을 곱한다는 것 = 선형 사상이라고 할 수 있음

 - 행렬끼리의 곱 = 사상의 합성

 - 영행렬$O$은 모든것을 영점으로 되돌리는 사상

 - 행렬의 거듭제곱 = 사상의 반복: $AA=A^{2}$과 같이 $A$가 공간에서의 사상이라고하면 2번 반복하는 의미가있다.

 - 행렬의 0제곱: $A^{0}=I$  I은 단위행렬을 의미

3. 단위행렬: 정방행렬에서 대각행렬만 1이고 모두 0인 행렬 $\begin{bmatrix}
 1 & 0 & 0\\ 
 0 & 1 & 0\\ 
 0 & 0 & 1
\end{bmatrix}$

4. 대각행렬: 정방행렬이지만, 대각행렬만 0이 아니고 모두 0인경우 $\begin{bmatrix}
 3 & 0 & 0\\ 
 0 & 1 & 0\\ 
 0 & 0 & 5
\end{bmatrix}$

5. 역행렬= 역사상. 정의: 정방행렬$A$에 대해 그 역사상에 대응하는 행렬을 A의 역행렬이라고 함. 기로호는 $A^{-1}$으로 쓰고, 어떠한 $x$을 가져와도 $Ax=y$ 또는 $A^{-1}y=x$이다. 그리고, 이 역행렬은 있을 수도 없을 수도 있음.

6. 행렬의 곱셈 트릭

$y=A\textbf{x}+b$ 을 행렬의 곱으로만 표시할 수 있을까? 정답은 Y이다. +b와 같이 상수항의 덧셈이 있을 땐, $\tilde{y} = \tilde{A}\tilde{x}$ 으로 표기하고, $\tilde{A}=A+o^{T}$, $\tilde{x}=x+1$ 으로 표기하면 한번에 처리할 수 있다.

7. 행렬식: determinant. 표기법: $det A, |A|$

8. 행렬식의 성질

 - det(I) = 1, 

 - det(AB) = det(A) x det(B)

 - $det A = \frac{1}{detA}$

 - 어느열이든 두 열이 완전히 같은 경우는 det A = 0. 예) $\begin{bmatrix}
 2 &  2& 4 \\ 
 7 &  7& 5 \\ 
 3 &  3& 2
\end{bmatrix}$

반응형

역전파알고리즘은 1986년 데이비드 등이 개발한 알고리즘이다. 핵심은 "은닉 유닛(hidden unit)이 무엇을 해야하는지는 모르지만, 은닉층의 활성도를 바꿀때 얼마나 빨리 오차가 변하는지는 계산할 수 있다".

 

그림1. Simple neural network

 

반환되는 Outcome $\hat{y}$가 변화할 때, 오차율(E)가 얼마나 변화할 것인지를 다음과 같이 미분식을 이용해보자. 아랫첨자 j은 데이터포인트이다. 

$\frac{\delta E}{\delta W_{1}} = \frac{\delta E}{\delta o_{1}} \frac{\delta o_{1}}{\delta z_{1}} \frac{\delta z_{1}}{\delta_{W_{1}}}$

 

첫번째 곱인 $\frac{\delta E}{\delta o_{1}}$ 은 MSE와 구하는 공식이 같아서 구할 수 있고, 

두번쨰 곱인 시그모이드(sigmoid)함수의 미분은 $\sigma'(\cdot) = \sigma(\cdot)(1-\sigma(\cdot))$ 으로 구할 수 있다.

세번째 곱은 W와 H1의 곱이 Z1이기 때문에,  $\frac{\delta z_{1}}{\delta_{W_{1}}} = h_{1}$

 

이를 모두 곱해주면된다. 곱한 모든 식의 아래와 같이 업데이트하는 과정을 이어나가면 끝이다.

$W_{1}^{new} = W_{1} -\alpha \frac{\delta E}{\delta W_{1}}$

 

 

위와 같은 과정을 $\hat{y}$ 을 이전레이어의 출력값으로 변경해서 층마다 업데이트하면된다.

 

 

예제를 위해서, 아래의 간단한 MLP(Multi-layer perceptron)을 생성해보았다.

import tensorflow as tf
model = tf.keras.models.Sequential()
model.add(tf.keras.layers.InputLayer(input_shape=(100)))
model.add(tf.keras.layers.Dense(10))
model.add(tf.keras.layers.Dense(1, activation='sigmoid'))

model.compile(loss='binary_crossentropy')

model.summary()

 

그리고 랜덤으로 xs, ys을 만들어주고, 학습시켰다.

xs = tf.constant(tf.random.uniform(shape=(100, 10)))
ys = tf.round(tf.random.uniform(shape=(100, 1)))

model.fit(xs, ys)

 

다행히 tensorflow2.x 은 역전파알고리즘을 tf.Gradient 객체로 한번에 구할 수 있게 만들어주어, 이를 사용하면된다.

with tf.GradientTape(persistent=True) as tape:
    tape.watch(xs)
    w = model.trainable_weights
    ys = model(xs)
    
    
dE_dW = tape.gradient(ys, w)

 

아래의 gradient descent방법과 같이 update을 해주면된다.

lr = 0.005
Ws = model.trainable_weights

for W, grad in zip(Ws, dE_dW):
    W = W - (lr * grad)

 

 

추가적으로 matrix derivation을 직접해보는경우 아래와같이 matrix shape이 안맞는 경우가 발생한다. 필자도 이 떄문에 한 일주일을 낭비했다..

 

마찬가지로, 우리가 구하고싶은것은 dE/dW 이므로, 아래와 같다.

 

 

https://medium.com/jun-devpblog/dl-3-backpropagation-98206058d72e

 

[DL] 3. Backpropagation

1. Error Backpropagation

medium.com

https://souryadey.github.io/teaching/material/Matrix_Calculus.pdf

반응형
import math
import numpy as np
import cv2

def ssim(img1, img2):
    C1 = (0.01 * 255)**2
    C2 = (0.03 * 255)**2

    img1 = img1.astype(np.float64)
    img2 = img2.astype(np.float64)
    kernel = cv2.getGaussianKernel(11, 1.5)
    window = np.outer(kernel, kernel.transpose())

    mu1 = cv2.filter2D(img1, -1, window)[5:-5, 5:-5]  # valid
    mu2 = cv2.filter2D(img2, -1, window)[5:-5, 5:-5]
    mu1_sq = mu1**2
    mu2_sq = mu2**2
    mu1_mu2 = mu1 * mu2
    sigma1_sq = cv2.filter2D(img1**2, -1, window)[5:-5, 5:-5] - mu1_sq
    sigma2_sq = cv2.filter2D(img2**2, -1, window)[5:-5, 5:-5] - mu2_sq
    sigma12 = cv2.filter2D(img1 * img2, -1, window)[5:-5, 5:-5] - mu1_mu2

    ssim_map = ((2 * mu1_mu2 + C1) * (2 * sigma12 + C2)) / ((mu1_sq + mu2_sq + C1) *
                                                            (sigma1_sq + sigma2_sq + C2))
    return ssim_map.mean()


def calculate_ssim(img1, img2):
    '''calculate SSIM
    the same outputs as MATLAB's
    img1, img2: [0, 255]
    '''
    if not img1.shape == img2.shape:
        raise ValueError('Input images must have the same dimensions.')
    if img1.ndim == 2:
        return ssim(img1, img2)
    elif img1.ndim == 3:
        if img1.shape[2] == 3:
            ssims = []
            for i in range(3):
                ssims.append(ssim(img1, img2))
            return np.array(ssims).mean()
        elif img1.shape[2] == 1:
            return ssim(np.squeeze(img1), np.squeeze(img2))
    else:
        raise ValueError('Wrong input image dimensions.')

 

원본(img1) 가우시안 필터(window), 적용한 이미지(mu1)

반응형

사용이유: 훈련시에 internal covariate shift

 

Batch noralization (BN)

 

훈련모드: 레이어가 미니배치 단위로 mean, std을 현재 미니배치단위마다 얻어와 각 채널별로 정규화해준다.

추론모드: 추론모드는 2가지 방법으로 사용될 수 있다. 

  • model.evaluate() 또는 prdict
  • layer, model이 training=False인 경우

추론모드에서는 각 훈련모드에서 얻어진 평균과 std을 moving averrage하여 번환한다.

  • moving_mean = moving_mean * momentum + mean(batch) * (1 - momentum)
  • moving_var = moving_var * momentum + var(batch) * (1 - momentum)

ma * (batch - self.moving_mean) / sqrt(self.moving_var + epsilon) + beta.

따라서, 추론모드는 이미 훈련모드에서 훈련된 상황에서만 사용할 수 있고, 훈련데이터가 추론데이터가 유사한 상황에서 사용될 수 있다. 따라서, 훈련모드와 추론모드에서 둘다 BN(Batch normalization)이 시행된다.

 

 

self.moving_mean_initializer = initializers.get(moving_mean_initializer)
self.moving_variance_initializer = initializers.get(moving_variance_initializer)


self.moving_mean = self.add_weight(
          name='moving_mean',
          shape=param_shape,
          dtype=self._param_dtype,
          initializer=self.moving_mean_initializer,
          synchronization=tf_variables.VariableSynchronization.ON_READ,
          trainable=False,
          aggregation=tf_variables.VariableAggregation.MEAN,
          experimental_autocast=False)


self.moving_variance = self.add_weight(
          name='moving_variance',
          shape=param_shape,
          dtype=self._param_dtype,
          initializer=self.moving_variance_initializer,
          synchronization=tf_variables.VariableSynchronization.ON_READ,
          trainable=False,
          aggregation=tf_variables.VariableAggregation.MEAN,
          experimental_autocast=False)
          
          
output, mean, variance = control_flow_util.smart_cond(training, train_op, _fused_batch_norm_inference)
 def mean_update():
   """Update self.moving_mean with the most recent data point."""
   
   if use_fused_avg_updates:
  		return self._assign_new_value(self.moving_mean, mean)
   else:
   		return self._assign_moving_average(self.moving_mean, mean, momentum,
   input_batch_size)

 

reference: https://keras.io/api/layers/normalization_layers/batch_normalization/

반응형
분석방법 적용 특징 문제 수식
Wilcox rank test        
KM method 각 군의자료가 적은경우   Time to event (발생할 시간)을 예측  
Cox proportional 제 3의 교란변수를 검정할 때 생존시간에 대한 별도의 가정이 없음+
공변량, 각 변수들이 주어졌을 때 식으로 표현이 가능함
= Semi-parametric 
Time to event (발생할 시간)을 예측 $\lambda(t|X_{i}) = \lambda_{0}(t)exp(X_{i}\beta)$

 

 

Cox proportional hazard model


각 i번째 환자(데이터포인트) X가 있다고 하자. $X_{i}=(X_{i1}, ...X_{ip})$ 와 같이 표기할 수 있으며 i번쨰 환자의 각 p라는 변수(공변량, covariates)를 의미한다. Cox비례위험모형은 아래와 같이 작성할 수 있다[식(2)]. 

 

$\lambda(t|X_{i}) = \lambda_{0}(t)exp(\beta_{1}X_{i1}+...+\beta_{p}X_{ip}$ [식2]

 

이를 Matrix 표기($X_{i}$)로 변경하면 흔히 보는 식이다[식3]. 보통 이를 hazard function이라고 한다.

 

 

 

$\lambda(t|X_{i}) = \lambda_{0}(t)exp(X_{i}\beta)$ [식3]

 

위의 식을 들여다보면, X가 커지든 작아지든 t보다는 $\beta$와 연관이 있다.

1) 즉, 시점이 고정이라면(환자들끼리 같은 조건이라면), X가 커질수록 $\beta$가 커져서, 위험이 올라간다고 할 수 있다. 즉 시간에 따라서는 위험에 대한 비율이 일정하게 유지가 된다(비례위험ㄱ자ㅓㅇ)

2) 또한, t에 영향을 받는것은 $\lambda_{0}$만 영향을 받기때문에, 시간이 증가할수록 커진다.

 

 

그리고, 흔히 위험비(harzard ratio, HR)라고 불리는 위험에 대한 비율(집단A에 대한 집단 B에 대한 이벤트 발생에 대한 비율)을 다음과 같이 구할 수 있다. 위의 hazard function을 A집단에 있는 환자a와 B집단에 있는 환자b에 대해서 적용하고자면 아래와 같이 각각 구할 수 있다.

 

$\lambda(t|X_{a}) = \lambda_{0}(t)exp(\beta_{1}X_{i1a}+...+\beta_{p}X_{ia})$

$\lambda(t|X_{b}) = \lambda_{0}(t)exp(\beta_{1}X_{i1a}+...+\beta_{p}X_{ib})$

 

위의 두 식을 나누면 아래와 같은 식을 구할 수 있는데, 보면 $\lambda(t)$은 없어지고, $X, \beta$로만 이뤄져있다. 즉, 시점 T에는 무관하게 위험이 동일하다고 할 수 있다. 

 

$\frac{\lambda(t|X_{a})}{\lambda(t|X_{b})} 
= \frac{exp(\beta_{1}X_{i1a}+...+\beta_{p}X_{ia})}{exp(\beta_{1}X_{i1a}+...+\beta_{p}X_{ib})} 
= exp(\beta(X_{a}-X_{b})) $

 

$HR=exp(\beta(X_{a}-X_{b}) = exp(\beta_{1}(X_{a1}-X_{b1})+...+beta_{p}(X_{ap}-X_{bp}))$을 의미하고, 각 $\beta_{p}$은 $X_{p}$에 따른 위험을 의미한다. 만약 다른 환자 a와 환자b가 모든 변수가 동일하고 $x_{1}$만 다르면 어떻게될까?

$\frac{\lambda(t|X_{a})}{\lambda(t|X_{b})} 

= \frac{h(t|x_{1}=1)}{h(t|x_{1}=0)}
= exp(\beta_{1}(1-0))
= exp(\beta_{1})
$

 

위처럼 모든 $x_{1}$ 변수만 제외하고 모든 변수가 동일하였을때 계수 $\beta_{1}$이 $x_{1}$가 1단위 변화할때마다 $exp(\beta)$한다고 해석하고, 이를 harzard ratio 라고 한다.

 

 

Python3에서는 lifelines패키지가 생존분석을 지원하는데, 시각화가 R보다는 약하다. 아래의 예시는 Cox을 이용한 시각화 및 통계플롯인데, CI값이 안나와서 다소 아쉽지만. python3에서 지원하는게 어딘가 싶다. (https://buildmedia.readthedocs.org/media/pdf/lifelines/latest/lifelines.pdf)

 

from lifelines.datasets import load_rossi
from lifelines import CoxPHFitter

rossi = load_rossi()
cph = CoxPHFitter()
cph.fit(rossi, duration_col='week', event_col='arrest')
cph.print_summary()


cph.plot()

cph.plot_partial_effects_on_outcome(covariates='race', values=list(set(rossi.race)), cmap='coolwarm')

 

반응형

+ Recent posts