Microsoft MVP성태의 닷넷 이야기
글쓴 사람
정성태 (techsharer at outlook.com)
홈페이지
첨부 파일

(시리즈 글이 7개 있습니다.)
Math: 15. 그래프 그리기로 알아보는 뉴턴-랩슨(Newton-Raphson's method)법과 제곱근 구하기 - C#
; https://www.sysnet.pe.kr/2/0/10911

Math: 53. C# - 행렬식을 이용한 최소 자승법(LSM: Least Square Method)
; https://www.sysnet.pe.kr/2/0/11918

Math: 54. C# - 최소 자승법의 1차 함수에 대한 매개변수를 단순 for 문으로 구하는 방법
; https://www.sysnet.pe.kr/2/0/11919

Math: 55. C# - 다항식을 위한 최소 자승법(Least Squares Method)
; https://www.sysnet.pe.kr/2/0/11921

Math: 56. C# - 그래프 그리기로 알아보는 경사 하강법의 최소/최댓값 구하기
; https://www.sysnet.pe.kr/2/0/11923

Math: 57. C# - 해석학적 방법을 이용한 최소 자승법
; https://www.sysnet.pe.kr/2/0/11924

Math: 58. C# - 최소 자승법의 1차, 2차 수렴 그래프 변화 확인
; https://www.sysnet.pe.kr/2/0/11936




C# - 최소 자승법의 1차 함수에 대한 매개변수를 단순 for 문으로 구하는 방법

일단 행렬식을 이용하면,

C# - 행렬식을 이용한 최소 자승법(LSM: Least Square Method)
; https://www.sysnet.pe.kr/2/0/11918

범용적으로 다항식에 대한 근사를 최소 자승법(최소 제곱법)으로 구할 수 있습니다. 하지만, 만약 대상을 "1차 함수"로 직선에 대한 근사만을 구한다면 복잡한 행렬 연산 없이 for 문만으로 매개 변수를 구하는 것이 가능합니다.

가령, 지난번 예제의 행렬식을 보겠습니다.

θ0 + θ1x1 = y1
θ0 + θ1x2 = y2
...
θ0 + θ1xn = yn



AX=B
A-1AX=A-1B
X=A-1B (A-1 == 의사역행렬)

결국 중요한 것은, 위의 식에서 A-1B 연산 결과를 구하는 것인데요, 헷갈리니까 일단 의사역행렬을 A+라고 정의하고, 이것을 행렬 라이브러리를 이용하면 단순히 Matrix 타입의 PseudoInverse를 호출하는 것으로 쉽게 해결했지만 만약 직접 구하고 싶다면 다음과 같은 과정을 거쳐야 합니다.

A+ = (ATA)-1AT

따라서, 매개변수를 나타내는 행렬 X는 B 행렬까지 곱해주면서 다음과 같이 계산할 수 있습니다.

X = A+ * B
  = (ATA)-1AT * B

이제 남은 작업은 위의 식을 간략하게 바꿔주면 됩니다. ^^




이 상태에서 A 행렬을 보면 "n x 2" 행렬이고 이것의 전치 행렬(AT)은 "2 x n" 행렬이 됩니다. 또한 B 행렬도 "n x 1" 행렬임을 감안하면 연산 결과가 다음과 같이 정리될 수 있습니다.

X = (ATA)-1AT * B
  = ((2 x n) * (n x 2))-1 * (2 x n) * (n x 1)
  = (2 x 2)-1 * (2 x 1)
  = (2 x 2) * (2 x 1)
  = (2 x 1)

즉, X 행렬은 (당연히 1차 함수의 매개변수 2개를 구하는 것이므로) 언제나 "2 x 1" 행렬이 나오므로 X 행렬의 인덱스에 해당하는 값을 정리해 볼 수도 있습니다. 이 과정을 단계별로 천천히 ^^ 접근해 볼까요?





2 x 2 행렬의 역행렬은 다음과 같이 간략화할 수 있으므로,



ATA 결과의 역행렬을 구할 수 있습니다.



위의 결과를 2 x n 행렬의 AT와 연산을 하면 과정이 좀 복잡하니 어차피 행렬곱은 결합법칙이 성립하므로 뒤의 AT B 연산을 먼저 다음과 같이 정리할 수 있습니다.



마지막으로 (ATA)-1(2 x 2 행렬)에 AT * B(2 x 1 행렬)을 곱하는 것이므로 다음과 같이 최종 정리가 됩니다.



따라서 1차 방정식의 매개변수는 이렇게 단일 식으로 각각 구할 수 있습니다.






구하는 과정에 정리할 식이 좀 끼어들어서 그렇지, 사실 C# 코드로 위의 계산을 나타내면 별거 아닙니다. ^^

// 단순 for 루프를 이용한 계산
private static (double theta1, double theta0) GetEquation2(double[] xData, double[] yData)
{
    double sumAnBn = 0.0;
    double sumAn = 0.0;
    double sumBn = 0.0;
    double sumAnAn = 0.0;

    for (int i = 0; i < xData.Length; i ++)
    {
        sumAnBn += xData[i] * yData[i];
        sumAn += xData[i];
        sumBn += yData[i];
        sumAnAn += xData[i] * xData[i];
    }

    int n = xData.Length;
    double Q = sumAnAn * n - sumAn * sumAn;

    double theta1 = (n * sumAnBn - sumAn * sumBn) / Q;
    double theta0 = (-sumAn * sumAnBn + sumAnAn * sumBn) / Q;

    return (theta1, theta0);
}

// 행렬을 이용한 계산
private static (double theta1, double theta0) GetEquation(double[] xData, double[] yData)
{
    Matrix matA = CreateMatrix.DenseOfColumnMajor(xData.Count(), 1, xData);
    Vector add1 = Vector.Build.DenseOfArray(Enumerable.Repeat(1.0, xData.Count()).ToArray());
    Matrix matAwith1 = matA.InsertColumn(1, add1);

    Console.WriteLine(matAwith1);
    Matrix matB = CreateMatrix.DenseOfColumnMajor(yData.Count(), 1, yData);

    Matrix pinvMatA = matAwith1.PseudoInverse();
    Console.WriteLine(pinvMatA);

    Matrix matX = pinvMatA * matB;

    return (matX[0, 0], matX[1, 0]);
}

당연하겠지만 행렬식을 이용했던 GetEquation 메서드와 비교해 보면,

{
    // y = theta0 + (theta1 * x)
    (double theta1, double theta0) = GetEquation(xData, yData);
    Console.WriteLine($"[method1] y = {theta0} + {theta1} * x");
    /* 출력 결과
    [method1] y = 231.545758451005 + 1.39551018043075 * x
    */
}

{
    (double theta1, double theta0) = GetEquation2(xData, yData);
    Console.WriteLine($"[method2] y = {theta0} + {theta1} * x");
    /* 출력 결과
    [method2] y = 231.545758451006 + 1.39551018043075 * x
    */
}

부동소수점 계산임을 감안해 값이 거의 동일하다는 것을 알 수 있습니다.

(첨부 파일은 이 글의 예제 코드를 포함합니다.)




[이 글에 대해서 여러분들과 의견을 공유하고 싶습니다. 틀리거나 미흡한 부분 또는 의문 사항이 있으시면 언제든 댓글 남겨주십시오.]







[최초 등록일: ]
[최종 수정일: 5/28/2019]

Creative Commons License
이 저작물은 크리에이티브 커먼즈 코리아 저작자표시-비영리-변경금지 2.0 대한민국 라이센스에 따라 이용하실 수 있습니다.
by SeongTae Jeong, mailto:techsharer at outlook.com

비밀번호

댓글 작성자
 



2019-05-28 10시35분
선형 최소 제곱법(Linear Least Squares Method) 사용하기
; https://icodebroker.tistory.com/5579
정성태

... 91  92  93  [94]  95  96  97  98  99  100  101  102  103  104  105  ...
NoWriterDateCnt.TitleFile(s)
11619정성태7/21/201823187Graphics: 10. Unity로 실습하는 Shader (6) - Mosaic Shading
11618정성태7/20/201819863개발 환경 구성: 387. 삼성 오디세이(Odyssey) 노트북의 운영체제를 새로 설치하는 방법
11617정성태7/20/201820487Team Foundation Server: 50. TFS 소스 코드 관리 기능 (5) - "Rollback", "Rollback Entire Changeset"
11616정성태7/17/201819569Graphics: 9. Unity Shader - 전역 변수의 초기화
11615정성태7/17/201823977.NET Framework: 788. RawInput을 이용한 키보드/마우스 입력 모니터링파일 다운로드1
11614정성태7/17/201826909Graphics: 8. Unity Shader - Texture의 UV 좌표에 대응하는 Pixel 좌표
11613정성태7/16/201823416Graphics: 7. Unity로 실습하는 Shader (5) - Flat Shading
11612정성태7/16/201821381Windows: 148. Windows - Raw Input의 Top level collection 의미
11611정성태7/15/201821748Graphics: 6. Unity로 실습하는 Shader (4) - 퐁 셰이딩(phong shading)
11610정성태7/15/201819147Graphics: 5. Unity로 실습하는 Shader (3) - 고로 셰이딩(gouraud shading) + 퐁 모델(Phong model) + Texture
11609정성태7/15/201822220Graphics: 4. Unity로 실습하는 Shader (2) - 고로 셰이딩(gouraud shading) + 퐁 모델(Phong model)
11608정성태7/15/201825845Graphics: 3. Unity로 실습하는 Shader (1) - 컬러 반전 및 상하/좌우 뒤집기
11607정성태7/14/201826218Graphics: 2. Unity로 실습하는 Shader [1]
11606정성태7/13/201827043사물인터넷: 19. PC에 연결해 동작하는 자신만의 USB 장치 만들어 보기파일 다운로드1
11605정성태7/13/201823398사물인터넷: 18. New NodeMCU v3 아두이노 호환 보드의 내장 LED 및 입력 핀 사용법 [1]파일 다운로드1
11604정성태7/12/201822396Math: 47. GeoGebra 기하 (24) - 정다각형파일 다운로드1
11603정성태7/12/201817641Math: 46. GeoGebra 기하 (23) - sqrt(n) 제곱근파일 다운로드1
11602정성태7/11/201818229Math: 45. GeoGebra 기하 (22) - 반전기하학의 원에 관한 반사변환파일 다운로드1
11601정성태7/11/201821259Math: 44. GeoGebra 기하 (21) - 반전기하학의 직선 및 원에 관한 반사변환파일 다운로드1
11600정성태7/10/201819914Math: 43. GeoGebra 기하 (20) - 세 점을 지나는 원파일 다운로드1
11599정성태7/10/201818825Math: 42. GeoGebra 기하 (19) - 두 원의 안과 밖으로 접하는 직선파일 다운로드1
11598정성태7/10/201821013Windows: 147. 시스템 복구 디스크를 USB 디스크에 만드는 방법
11597정성태7/10/201823252사물인터넷: 17. Thinary Electronic - ATmega328PB 아두이노 호환 보드의 개발 환경 구성
11596정성태7/10/201820647기타: 72. 과거의 용어 설명 - OWIN
11595정성태7/10/201826140사물인터넷: 16. New NodeMCU v3 아두이노 호환 보드의 기본 개발 환경 구성
11594정성태7/8/201821085Math: 41. GeoGebra 기하 (18) - 원의 중심 및 접선파일 다운로드1
... 91  92  93  [94]  95  96  97  98  99  100  101  102  103  104  105  ...