데이터마이너를 꿈꾸며

보스턴 하우징을 이용한 신경망 모형. 목표변수를 주택가격의 중간값인 medv(연속형 변수)로 하고, 나머지 변수를 입력변수로 하는 신경망 모형

1) 데이터 입력

> set.seed(100)

> library(MASS)

> library(neuralnet)

> bdat = Boston

> str(bdat)

'data.frame':     506 obs. of  15 variables:

$ crim    : num  0.00632 0.02731 0.02729 0.03237 0.06905 ...

$ zn      : num  18 0 0 0 0 0 12.5 12.5 12.5 12.5 ...

$ indus   : num  2.31 7.07 7.07 2.18 2.18 2.18 7.87 7.87 7.87 7.87 ...

$ chas    : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...

$ nox     : num  0.538 0.469 0.469 0.458 0.458 0.458 0.524 0.524 0.524 0.524 ...

$ rm      : num  6.58 6.42 7.18 7 7.15 ...

$ age     : num  65.2 78.9 61.1 45.8 54.2 58.7 66.6 96.1 100 85.9 ...

$ dis     : num  4.09 4.97 4.97 6.06 6.06 ...

$ rad     : Factor w/ 9 levels "1","2","3","4",..: 1 2 2 3 3 3 5 5 5 5 ...

$ tax     : num  296 242 242 222 222 222 311 311 311 311 ...

$ ptratio : num  15.3 17.8 17.8 18.7 18.7 18.7 15.2 15.2 15.2 15.2 ...

$ black   : num  397 397 393 395 397 ...

$ lstat   : num  4.98 9.14 4.03 2.94 5.33 ...

$ medv    : num  24 21.6 34.7 33.4 36.2 28.7 22.9 27.1 16.5 18.9 ...

$ medv.hat: num  23.7 23.7 32.4 32.4 32.4 ...

 

주의: chas, rad를 수치형 변수로 바꾸고, 앙상블 모형에서 사용하였던 medv.hat 변수를 삭제해야 한다.

신경망에서는 범주형 데이터를 수치화하여 적용한다. 대체로 순서가 의미 있는 범부형 데이터는 수치 전환을 한 뒤 표준화하여 사용하게 된다.

 

2) 데이터 및 타입 변경

> bdat<-bdat[,-15]

> bdat$chas = as.numeric(bdat$chas)

> bdat$rad = as.numeric(bdat$rad)

> class(bdat$chas);class(bdat$rad)

[1] "numeric"

[1] "numeric"

 

3) 50% 랜덤 추출

> i = sample(1:nrow(bdat), round(0.5*nrow(bdat)))

> max1 = apply(bdat, 2, max)

> min1 = apply(bdat, 2, min)

결과 해석 : 50% 랜덤 추출하여 훈련데이터(train)로 저장하고 나머지는 검증 데이터(test)로 저장

 

 

4) 변수의 표준화 과정

> sdat = scale(bdat,center=min1,scale = max1 - min1) 

> sdat = as.data.frame(sdat)

> head(sdat,3)

crim          zn               indus  chas                nox                   rm                 age                   dis   rad

1 0.0000000000000 0.18 0.06781524927    0 0.3148148148 0.5775052692 0.6416065911 0.2692031391 0.000

2 0.0002359225392 0.00 0.24230205279    0 0.1728395062 0.5479977007 0.7826982492 0.3489619802 0.125

3 0.0002356977440 0.00 0.24230205279    0 0.1728395062 0.6943858977 0.5993820803 0.3489619802 0.125

tax               ptratio              black                     lstat         medv

1 0.2080152672 0.2872340426 1.0000000000 0.08967991170 0.4222222222

2 0.1049618321 0.5531914894 1.0000000000 0.20447019868 0.3688888889

3 0.1049618321 0.5531914894 0.9897372535 0.06346578366 0.6600000000

 

함수 설명:

> sdat = scale(bdat,center=min1,scale = max1 - min1) 신경망에서 사용하는 수치형 변수를 0과 1사이의 값으로 바꾸어 주기 위한 단계. 위의 [0-1] 변환은 (x-min(x) /(max(x)-min(x)) 의 수식으로도 계산 가능. 

 

 

5) 신경망 모델 구축 및 신경망 그래프 그리기

> train = sdat[i,] #학습, 훈련샘플 training

> test = sdat[-i,] #테스트샘플 test

> n = names(train) ;n

[1] "crim"    "zn"      "indus"   "chas"    "nox"     "rm"      "age"     "dis"     "rad"     "tax"   

[11] "ptratio" "black"   "lstat"   "medv"  

> form = as.formula(paste('medv~',paste(n[!n %in% 'medv'],collapse = '+'))) 

> form

medv ~ crim + zn + indus + chas + nox + rm + age + dis + rad + tax + ptratio + black + lstat

> nn1 = neuralnet(form,data=train,hidden = c(5,3),linear.output = T)

> summary(nn1) #작성된 신경망모형의 오브젝트의 구조를 정리

                       Length    Class     Mode   

call                          5   -none-     call   

response              253   -none-     numeric

covariate            3289   -none-     numeric

model.list                 2   -none-     list   

err.fct                      1   -none-     function

act.fct                      1   -none-     function

linear.output             1   -none-     logical

data                       14   data.frame list   

net.result                 1   -none-     list   

weights                    1   -none-     list   

startweights              1   -none-     list   

generalized.weights   1   -none-     list   

result.matrix            95   -none-     numeric

함수 설명

> form = as.formula(paste('medv~',paste(n[!n %in% 'medv'],collapse = '+'))) #종속변수는 mdev, 나머지는 독립변수이다 틸다 ~. 와같은 개념. 종속변수를 제외한 n 안에 있는 모든 변수명을 다 집어넣고 더하여라. 변수명을 직접 입력해도 된다.

실제 form라고 입력하면 medv ~ crim + zn + indus + chas + nox + rm + age + dis + rad + tax + ptratio + black + lstat

 

라는 공식을 볼 수 있다.

 

> nn1 = neuralnet(form,data=train,hidden = c(5,3),linear.output = T) #은닉층은 2개이고 처음은 5개 나머지는 3개, 회귀의 문제이므로 linear.output = T

 

 

질문: 은닉층의 수는 어떻게 결정할까? 통상 은닉층의 마디수는 입력층의 마디수 두배를 넘지 않도록 해야한다. 위에서는 c(5,3)으로 결정하였는데 이 수치가 최적의 수치인지는 어떻게 알까? 더 공부가 필요하다.

 

<참고> 딥러닝(deep learning)

기계학습 기법 중 하나로 신경망모형으로부터 비롯된 딥러닝(deep learning)은 기본적으로 은닉층이 많이 쌓여 가면서 복잡하고 깊은 구조로 발전하면서 deep 이라는 이름이 붙여졌다. 입력변수와 출력변수 간 복잡한 관계를 가중치를 통해 조정할 수 있는 구조와 매커니즘을 갖고 있다.

 

 

> plot(nn1)

 

 

①로 표시된 것이 상수항에 해당하고 가중치는 각각의 화살 표 위에 출력된다. 해석 공부도 더 필요하다.

처음 은닉층은 5개, 두번째 은닉층은 3개를 가지고 있는 신경망 모형이 완성된다. 

6) 모형 추정: 위의 그래프는 50%의 training data만을 사용하였으니 나머지 50% test data를 쓰자.

> pred.nn0 = compute(nn1,test[,1:13])

> summary(pred.nn0)
               Length Class  Mode  
    neurons      3    -none- list  
    net.result 253    -none- numeric

> names(pred.nn0)

[1] "neurons"    "net.result"

> pred0 = pred.nn0$net.result*(max(bdat$medv)-min(bdat$medv))+min(bdat$medv)

> head(pred0)
                    [,1]    
    1 26.55282080
    2 25.15912380
    3 34.48852794
    4 31.56317828
    5 32.47594868
    6 25.32302689

 

함수 설명:

> pred.nn0 = compute(nn1,test[,1:13])

14번째 변수가 medv 목표변수.compute는 작성된 신경망모형을 이용하여 새로운 예에 적용하여 결과를 도출. nn1는 새롭게 예측에 적용할 자료, test[,1:13]는 신경망모형으로 적합한 결과 오브젝트

> pred0 = pred.nn0$net.result*(max(bdat$medv)-min(bdat$medv))+min(bdat$medv)

목표변수를 조정전 변수값으로 환원하는 과정. 원래 가지고 있는 값으로 환원. 역변환 과정

 

 

7) 학습샘플의 실제값과 예측값 (적합값)을 비교해보자. 

아래 함수 head(cbind(bdat[i,14],pred0),50)에서 cbind(bdat[i,14],는 실제값(training data)과 예측값(neural network)으로 구성된 행렬을 구성하는 함수

> head(cbind(bdat[i,14],pred0),50)

[,1]        [,2]

1  15.6 28.12482177

2  19.2 21.36405177

3  29.1 34.24035993

4  18.4 35.88436410

6  24.0 25.77273199

9  22.2 14.86779840

12 11.9 20.82595913

14 26.4 20.03825755

18 24.4 17.19077183

19 23.9 18.69934412

21 20.3 14.85391889

23  9.6 15.77921679

24 13.3 14.95074253

25 33.3 16.19467498

26 15.0 15.53094329

27 19.3 16.50978025

28 27.5 15.65043905

30 22.6 19.71434020

33 29.4 13.39051466

34 19.5 15.45310394

36 33.8 22.75202201

41 31.2 37.96410157

42 21.2 36.32787047

43 19.9 27.29495234

44 42.3 26.74797523

46 24.8 20.55343697

47 50.0 19.83522653

48 20.8 16.66204238

50 48.8 16.90901569

51 23.0 19.54979162

54 41.7 22.43724314

55 17.1 16.69040814

57 25.0 22.91956778

58 15.2 29.22726384

62  8.1 18.82030494

63  8.8 24.13500281

66 19.7 25.15181243

67 28.6 20.16650282

68 20.2 20.27218139

69 18.7 17.48036500

70 17.0 19.87776814

71 12.1 25.21432307

73 25.0 23.16314867

75 12.3 23.82067236

77 23.9 20.42068536

79 37.6 20.09872166

80 22.7 20.89416546

81  5.0 28.30078783

84 28.4 22.10685814

87 14.0 21.30514814

결과 해석: 26번째 행의 값처럼 유사한 결과값도 있지만 81번 결과 값처럼 서로 차이가 나는 경우도 있다.

 

 

8) 예측 정확도 평가: 모형 추정에 사용되지 않은 test data를 가지고 예측을 해보자. 

> pred.nn1 = compute(nn1,test[,1:13]) #목표 변수 제외, 새로운 변수로 간주하고 nn1에 적용

> pred1 = pred.nn1$net.result*(max(bdat$medv)-min(bdat$medv))+min(bdat$medv) #목표변수를 조정전 변수값으로 환원. 역변환

> head(cbind(bdat[-i,14],pred1),50) #좌측은 test 변수의 실제값, 우측은 예측값(적합값), 14번째 값은 목표변수, pred1은 예측값(적합값)

[,1]        [,2]

1  24.0 28.12482177

2  21.6 21.36405177

3  34.7 34.24035993

4  33.4 35.88436410

6  28.7 25.77273199

9  16.5 14.86779840

12 18.9 20.82595913

14 20.4 20.03825755

18 17.5 17.19077183

19 20.2 18.69934412

21 13.6 14.85391889

23 15.2 15.77921679

24 14.5 14.95074253

25 15.6 16.19467498

26 13.9 15.53094329

27 16.6 16.50978025

28 14.8 15.65043905

30 21.0 19.71434020

33 13.2 13.39051466

34 13.1 15.45310394

36 18.9 22.75202201

41 34.9 37.96410157

42 26.6 36.32787047

43 25.3 27.29495234

44 24.7 26.74797523

46 19.3 20.55343697

47 20.0 19.83522653

48 16.6 16.66204238

50 19.4 16.90901569

51 19.7 19.54979162

54 23.4 22.43724314

55 18.9 16.69040814

57 24.7 22.91956778

58 31.6 29.22726384

62 16.0 18.82030494

63 22.2 24.13500281

66 23.5 25.15181243

67 19.4 20.16650282

68 22.0 20.27218139

69 17.4 17.48036500

70 20.9 19.87776814

71 24.2 25.21432307

73 22.8 23.16314867

75 24.1 23.82067236

77 20.0 20.42068536

79 21.2 20.09872166

80 20.3 20.89416546

81 28.0 28.30078783

84 22.9 22.10685814

87 22.5 21.30514814

 

결과해석: 전반적으로 실제값과 예측값이 서로 유사함을 알 수 있다.

 

> head(cbind(bdat[-i,14],pred1)[,1],10)
   1    2    3    4    5    6    7    8    9   10
24.0 21.6 34.7 33.4 36.2 28.7 22.9 27.1 16.5 18.9
> head(cbind(bdat[-i,14],pred1),10)
   [,1]        [,2]
1  24.0 26.55282080
2  21.6 25.15912380
3  34.7 34.48852794
4  33.4 31.56317828
5  36.2 32.47594868
6  28.7 25.32302689
7  22.9 19.92609580
8  27.1 20.01203029
9  16.5 19.68915018
10 18.9 19.33552651
> obs <-cbind(bdat[-i,14],pred1)[,1]
> pdt <-cbind(bdat[-i,14],pred1)[,2]

> out = cbind(obs,pdt)
> head(out)
   obs         pdt
1 24.0 26.55282080
2 21.6 25.15912380
3 34.7 34.48852794
4 33.4 31.56317828
5 36.2 32.47594868
6 28.7 25.32302689
> which.min(abs(out$obs - out$pdt))
Error in out$obs : $ operator is invalid for atomic vectors
> out = as.data.frame(out)
> which.min(abs(out$obs - out$pdt))
[1] 43
> out[43,]
    obs         pdt
43 25.3 25.31048382
> which.max(abs(out$obs - out$pdt))
[1] 197
> out[197,]
    obs         pdt
372  50 18.13253008

 

 

9) PMSE: 예측된 값이 얼마나 잘 예측되었을까? 

> PMSE = sum((bdat[-i,14] - pred1)^2) / nrow(test)

> PMSE

[1] 18.70948041

함수 설명:

> PMSE = sum((bdat[-i,14] - pred1)^2) / nrow(test) #예측 평균 제곱 오차, 

 

 결과 해석: 예측평균제곱오차는 18.7

출처: 데이터마이닝(장영재, 김현중, 조형준 공저)