다중 선형 회귀를 iteration on data로 푸는 방법

 

자료

 

1111120
1111226
1111332
1112125
1112231
1112337
1113130
1113236
1113342
1121124
1121230
1121336
1122129
1122235
1122341
1123134
1123240
1123346
1131128
1131234
1131340
1132133
1132239
1132345
1133138
1133244
1133350
1211123
1211229
1211335
1212128
1212234
1212340
1213133
1213239
1213345
1221127
1221233
1221339
1222132
1222238
1222344
1223137
1223243
1223349
1231131
1231237
1231343
1232136
1232242
1232348
1233141
1233247
1233353
1311126
1311232
1311338
1312131
1312237
1312343
1313136
1313242
1313348
1321130
1321236
1321342
1322135
1322241
1322347
1323140
1323246
1323352
1331134
1331240
1331346
1332139
1332245
1332351
1333144
1333250
1333356
2111122
2111228
2111334
2112127
2112233
2112339
2113132
2113238
2113344
2121126
2121232
2121338
2122131
2122237
2122343
2123136
2123242
2123348
2131130
2131236
2131342
2132135
2132241
2132347
2133140
2133246
2133352
2211125
2211231
2211337
2212130
2212236
2212342
2213135
2213241
2213347
2221129
2221235
2221341
2222134
2222240
2222346
2223139
2223245
2223351
2231133
2231239
2231345
2232138
2232244
2232350
2233143
2233249
2233355
2311128
2311234
2311340
2312133
2312239
2312345
2313138
2313244
2313350
2321132
2321238
2321344
2322137
2322243
2322349
2323142
2323248
2323354
2331136
2331242
2331348
2332141
2332247
2332353
2333146
2333252
2333358
3111124
3111230
3111336
3112129
3112235
3112341
3113134
3113240
3113346
3121128
3121234
3121340
3122133
3122239
3122345
3123138
3123244
3123350
3131132
3131238
3131344
3132137
3132243
3132349
3133142
3133248
3133354
3211127
3211233
3211339
3212132
3212238
3212344
3213137
3213243
3213349
3221131
3221237
3221343
3222136
3222242
3222348
3223141
3223247
3223353
3231135
3231241
3231347
3232140
3232246
3232352
3233145
3233251
3233357
3311130
3311236
3311342
3312135
3312241
3312347
3313140
3313246
3313352
3321134
3321240
3321346
3322139
3322245
3322351
3323144
3323250
3323356
3331138
3331244
3331350
3332143
3332249
3332355
3333148
3333254
3333360

x0 ~ x5까지, 그리고 y

b0 ~ b5까지 추정해야 함

 

위 자료를 iod_multiple_linear_regression.dat로 저장

 

LHS와 RHS를 만드는 프로그램

 

PROGRAM iod_mlr_setup

! program name
! : iteration on data of multiple linear regreesion
! prgrammer
! : Park Byoungho
! date
! : 2009.11.1.

IMPLICIT NONE

! data dictionary

INTEGER, PARAMETER :: no_of_rc = 5 ! number of regression coefficient
CHARACTER(LEN=256), PARAMETER :: data_filename = 'iod_multiple_linear_regression.dat' ! data file name

REAL(KIND = 8), ALLOCATABLE :: ind(:) ! storage for independent variable when reading data
REAL(KIND = 8) :: dep ! storage for dependent variable when reading data

REAL(KIND = 8), ALLOCATABLE :: rhs(:) ! right-hand side
INTEGER, ALLOCATABLE :: eq_no(:) ! equation number

REAL(KIND = 8) :: temp_ind ! temporary independent variable
INTEGER :: temp_eq_no

INTEGER :: i,j ! for loop
INTEGER :: status ! I/O status
CHARACTER(LEN = 40) :: error_msg ! error message

ALLOCATE(ind(no_of_rc)) ! independent variable
ALLOCATE(rhs(no_of_rc)) ! right-hand side
ALLOCATE(eq_no(no_of_rc))

! initialize
rhs = 0

! open data file
OPEN(UNIT = 10, FILE = data_filename, STATUS = 'OLD', ACTION = 'READ', IOSTAT = status, IOMSG = error_msg)

IF (status /= 0) THEN ! file open failed
WRITE (*,'(1X, A, A)') 'Data file open failed -- error message : ', error_msg
STOP
END IF

! open data file for writing left-hand side
OPEN(UNIT = 20, FILE = 'lhs.dat', STATUS = 'REPLACE', ACTION = 'WRITE')

! read each line
DO
READ(UNIT = 10, FMT = *, IOSTAT = status) ind, dep
IF (status /= 0) EXIT ! end of file

DO i = 1, no_of_rc
eq_no(i) = i
END DO

DO i = 1, no_of_rc
! swap equation number
temp_eq_no = eq_no(i)
eq_no(i) = eq_no(1)
eq_no(1) = temp_eq_no

! swap independent variable
temp_ind = ind(i)
ind(i) = ind(1)
ind(1) = temp_ind

WRITE(20,*) (eq_no(j),ind(1)*ind(j), j=1,no_of_rc)

rhs(i) = rhs(i) + ind(1) * dep
END DO

END DO

CLOSE(11) ! close data file

! open file for right hand side
OPEN(UNIT = 30, FILE = 'rhs.dat', STATUS = 'REPLACE', ACTION = 'WRITE')

! write right-hand side
DO i = 1, no_of_rc
WRITE(30, *) rhs(i)
END DO

CLOSE(30) ! close

END PROGRAM iod_mlr_setup

 

위 프로그램 소스를 iod_multiple_linear_regression_setup.f95로 저장

 

LHS의 정렬

 

sort -n -o sorted_lhs.dat lhs.dat

 

정렬된 LHS와 RHS를 이용하여 해를 구하는 프로그램

 

PROGRAM iod_mlr_solve

! program name
! : iteration on data of multiple linear regreesion
! : solve program
! prgrammer
! : Park Byoungho
! date
! : 2009.11.1.

IMPLICIT NONE

! data dictionary

INTEGER, PARAMETER :: no_of_rc = 5 ! number of regression coefficient

INTEGER :: no_of_lhs ! number of left hand side lines
REAL(KIND = 8), ALLOCATABLE :: lhs(:,:), rhs(:) ! left-hand side, right-hand side, x'y
REAL(KIND = 8), ALLOCATABLE :: solutions(:) ! solutions

INTEGER :: pre_eq_no ! previous equation number
INTEGER :: cur_eq_no ! current equation number

REAL(KIND = 8) :: diag_ele ! diagonal element
REAL(KIND = 8) :: temp_rhs ! temporary rihgt hand side
REAL(KIND = 8) :: pre_sol ! previous solution

INTEGER :: iteration ! iteration
REAL(KIND = 8) :: epsilon ! sum of squares of difference between old and new solutions
INTEGER :: i,j ! loop
INTEGER :: status ! i/o status
CHARACTER(LEN = 40) :: error_msg ! error message
INTEGER, PARAMETER :: MAX_ITER = 500 ! maximum number of iteration
REAL(KIND = 8), PARAMETER :: criteria = 1.E-12 ! criteria for stopping

ALLOCATE(rhs(no_of_rc))
ALLOCATE(solutions(no_of_rc))

solutions = 0

!open left-hand side file
OPEN(UNIT = 10, FILE = 'sorted_lhs.dat', STATUS = 'OLD', ACTION = 'READ', IOSTAT = status, IOMSG = error_msg)
IF (status /= 0) THEN ! file open failed
WRITE (*,'(1X, A, A)') 'Sorted lhs file open failed -- error message : ', error_msg
STOP
END IF

! count the lines of lhs
no_of_lhs = 0
DO
READ(UNIT = 10, FMT = *, IOSTAT = status)
IF (status /= 0) EXIT ! reach end of file
no_of_lhs = no_of_lhs + 1
END DO

ALLOCATE(lhs(no_of_lhs, no_of_rc * 2))

! write(*,*) 'number of lines of left hand side = ', no_of_lhs

! store lhs to array
REWIND(10)
DO i = 1, no_of_lhs
READ(UNIT = 10, FMT = *, IOSTAT = status) lhs(i,:)
IF (status /= 0) EXIT ! reach end of file
END DO
CLOSE(11)

! DO i = 1, no_of_lhs
! write(*,*) lhs(i,:)
! END DO

!open right-hand side file
OPEN(UNIT = 20, FILE = 'rhs.dat', STATUS = 'OLD', ACTION = 'READ', IOSTAT = status, IOMSG = error_msg)

IF (status /= 0) THEN ! file open failed
WRITE (*,'(1X, A, A)') 'RHS file open failed -- error message : ', error_msg
STOP
END IF

! read and store right-hand side
DO i = 1, no_of_rc
READ(UNIT = 20,FMT = *) rhs(i)
END DO
close(20)

! write(*,*) 'right hand side ...'
! DO i = 1, no_of_rc
! write(*,*) i, rhs(i)
! END DO

DO iteration = 1, MAX_ITER

pre_eq_no = 1
diag_ele = 0.0
temp_rhs = rhs(1)
epsilon = 0.0

! start reading and processing each line of lhs
DO i = 1, no_of_lhs

cur_eq_no = lhs(i,1)

IF(pre_eq_no /= cur_eq_no) THEN
pre_sol = solutions(pre_eq_no) ! store old solution
solutions(pre_eq_no) = temp_rhs / diag_ele ! get a new solution
epsilon = epsilon + (pre_sol - solutions(pre_eq_no)) ** 2 ! calculate sum of squares of difference between old solution and new solution

diag_ele = 0.0
temp_rhs = rhs(INT(lhs(i,1)))
END IF

! add diagonal element
diag_ele = diag_ele + lhs(i, 2)

! adjust right-hand sie
! temp_rhs - independend variable * solutions
DO j = 2, no_of_rc
temp_rhs = temp_rhs - lhs(i, 2*j) * solutions(INT(lhs(i,2*j-1)))
END DO

pre_eq_no = lhs(i,1)

END DO
! end reading lhs

! calculate last solution
pre_sol = solutions(pre_eq_no) ! store old solution
solutions(pre_eq_no) = temp_rhs / diag_ele ! get a new solution
epsilon = epsilon + (pre_sol - solutions(pre_eq_no)) ** 2 ! calculate sum of squares of difference between old solution and new solution

! WRITE(*,*) iteration,'th iteration''s solutions'
! DO i = 1, no_of_rc
! WRITE(*,*) i, solutions(i)
! END DO

epsilon = epsilon / no_of_rc

! write iteration number and epsilon
write(*,*) 'iteration = ', iteration , ', epsilon = ', epsilon
IF (epsilon < criteria) THEN
EXIT
END IF

END DO
! end iteration

! open file for writing solution
OPEN(UNIT=30, FILE='sol.dat', STATUS='REPLACE', ACTION='WRITE', IOSTAT=status)
DO i = 1, no_of_rc
WRITE(30,*) i, solutions(i)
END DO
CLOSE(30)

END PROGRAM iod_mlr_solve

 

위 프로그램 소스를 iod_multiple_linear_regression_solve.f95로 저장

 

프로그램 컴파일 및 실행

 

 









관련 프로그램

 

1257082823_iod_multiple_linear_regression.zip
다운로드

1257082823_iod_multiple_linear_regression.zip

 

 

 

+ Recent posts