Excel 매크로 (Excel VBA)로 라그랑주 보간

소개

엑셀로 라그랑주 보간을 할 필요가 있었으므로, 매크로를 작성했습니다.
함수를 지정하고, 원하는 지점과 데이터를 지정하는 것만으로 간단하게 보간을 할 수 있습니다.
라그랑주 보간은 점수를 많게 하면 해의 진동이 발생하기 때문에, 표준에서는 근방 4점만으로 보간을 실시하도록 하고 있습니다.

함수 사양

함수 이름: LagrangeInterpolation(x, 알려진 y, 알려진 x, 이웃만 보간할 것인가)
인수:
x: 원하는 위치
알려진 y: 데이터 범위(Y)
알려진 x: 데이터 범위(X)
근방만: False로 하면 모든 지점에서 보간한다. 지정하지 않으면 자동으로 True

※인수가 "x""알려진 y""알려진 x"의 순서로 되어 있는 것은 FORECAST 함수와 같은 인수로 하고 싶기 때문입니다.

출처

Lagrange 보간.bas

Option Explicit

Function LagrangeInterpolation(X As Double, 既知のy As Range, 既知のx As Range, Optional 近傍のみ As Boolean = True) As Variant

    ' 近傍のみで補間する場合は前後N/2点ずつ補間する
    ' 例：N=4の場合は前後2点ずつ
    ' Nは2の倍数である必要がある
    Const N = 4

    ' yとxのデータは同じ個数であることが前提
    If 既知のx.count <> 既知のy.count Then
        ' xとyでデータ個数が異なる場合はN/Aを返す
        LagrangeInterpolation = CVErr(xlErrNA)
        Exit Function
    End If

    ' 生データの個数を代入
    Dim DataNum As Long
    DataNum = 既知のy.count

    ' ループ関数を定義
    Dim i As Long
    Dim j As Long

    ' Rangeデータを配列に保存
    Dim Xdatas() As Double
    ReDim Xdatas(既知のx.count)
    Dim Ydatas() As Double
    ReDim Ydatas(既知のy.count)
    For i = 1 To 既知のx.count
        Xdatas(i) = 既知のx(i)
        Ydatas(i) = 既知のy(i)
    Next

    ' 地点をXでソートし、X前後の地点をN点絞り込む
    ' 地点数Nよりデータ数のほうが少ない場合は何もしない
    If 近傍のみ And N <= DataNum Then

        ' ソートを行う
        Dim BufDouble As Double
        For i = 1 To UBound(Xdatas)
            For j = 1 To UBound(Xdatas)
                If Xdatas(i) < Xdatas(j) Then
                    BufDouble = Xdatas(i)
                    Xdatas(i) = Xdatas(j)
                    Xdatas(j) = BufDouble

                    BufDouble = Ydatas(i)
                    Ydatas(i) = Ydatas(j)
                    Ydatas(j) = BufDouble
                End If
            Next
        Next

        ' 切り出した後のXとYのデータ
        Dim CuttedXDatas() As Double
        ReDim CuttedXDatas(N)
        Dim CuttedYDatas() As Double
        ReDim CuttedYDatas(N)

        ' Xに最も近い要素番号を求める
        Dim StartCutNum As Long
        If X < Xdatas(2) Then
            ' Xより小さいデータが十分に存在しない場合、1～Nまでで補間を行う
            StartCutNum = 1
        ElseIf Xdatas(UBound(Xdatas) - 2) < X Then
            ' Xより大きいデータが十分に存在しない場合、生データ数-N～生データ数までで補間を行う
            StartCutNum = UBound(Xdatas) - N + 1
        Else
            ' X近傍に十分なデータ数がある場合、前後で補間
            For i = 1 To UBound(Xdatas) - 1
                If Xdatas(i) < X And X <= Xdatas(i + 1) Then
                    StartCutNum = i
                    Exit For
                End If
            Next
        End If

        If StartCutNum = 1 Then
            ' 1～N点で補間を行う
            For i = 1 To N
                CuttedXDatas(i) = Xdatas(i)
                CuttedYDatas(i) = Ydatas(i)
            Next
        ElseIf StartCutNum = UBound(Xdatas) - N + 1 Then
            ' 生データ数-N～生データ数で補間を行う
            For i = StartCutNum To DataNum
                CuttedXDatas(i - StartCutNum + 1) = Xdatas(i)
                CuttedYDatas(i - StartCutNum + 1) = Ydatas(i)
            Next
        Else
            ' 前後の地点で補間を行う
            For i = StartCutNum - N / 2 To StartCutNum + 1
                CuttedXDatas(i - StartCutNum + N / 2 + 1) = Xdatas(i)
                CuttedYDatas(i - StartCutNum + N / 2 + 1) = Ydatas(i)
            Next
        End If

        ' 切り出したデータを反映
        DataNum = N
        ReDim Xdatas(N)
        Xdatas = CuttedXDatas
        ReDim Ydatas(N)
        Ydatas = CuttedYDatas

    End If

    ' ラグランジュ既定関数の分子および分母
    Dim Numerator As Double     ' 分子
    Dim Denominator As Double   ' 分母

    ' 関数の戻り値を初期化
    LagrangeInterpolation = 0

    ' ラグランジュ既定関数を計算し、yの値との積を合計する
    For i = 1 To DataNum
        Numerator = 1
        Denominator = 1
        For j = 1 To DataNum
            If i <> j Then
                Numerator = Numerator * (X - Xdatas(j))
                Denominator = Denominator * (Xdatas(i) - Xdatas(j))
            End If
        Next
        LagrangeInterpolation = LagrangeInterpolation + (Ydatas(i) * Numerator / Denominator)
    Next

End Function

사용 예 및 결과

이 함수를 사용하여 다음 수식을 보간해 봅시다.

\frac{1}{1+8x^2}

엑셀 스크린샷

보간 후 그래프

검은 점선: 참값
얇은 붉은 선: 모든 점에서 보간된 결과(해가 진동함)
薄青太線:근방 4점에서 보간한 결과(해가 진동하지 않는다)

요약

데이터 수가 많으면 단순히 라그랑주 보간을 하는 것만으로는 해가 진동해 버립니다.
그 때문에 이번은 표준으로, 근방 지점만으로 보간을 하도록(듯이) 고안했습니다.
다만, 근방 지점만으로 보간해도 진동을 피할 수 없는 경우도 있으므로, 해의 취급에는 주의가 필요합니다.