係数に制約をつけてfittingする方法

Question

ExcelでもPythonでもRでもgnuplotでもまあツールはなんでも良いのですが、今3次式近似を行おうとしています。

図はExcelでやらした結果です。

ここで近似線は

y = ax^3+bx^2+cx+d

という式になっているわけですが、この式が(0,0), (1,1)を必ず通るように制約を書けたいのです。

つまり

d = 0

a + b + c = 1

という感じでしょうか。

制約なしでfittingする方法はいくらでも出てきますが、3次式近似に制約を掛ける例が見つけられませんでした。

上記のようなグラフと係数、R値を出す方法を教えて下さい。

Answer 1

直接的な回答ではありませんがコメントします。

被説明変数がモル分率なので、被説明変数は0～1の範囲であると制約を解釈し直します。この説明変数は任意の実数（実際には0～1の範囲内）で被説明変数は0～1の範囲の実数を満たす関数にLogit関数があります。そこでこの関数を応用してモデル化すると要望に近い結果となる結果が得られると考えます。

ただし、Logit関数は右肩上がりのS字曲線を描くのでグラフのデータをy=x(グラフ上の実線)に対して線対称に変換する必要があります。

以上をまとめると以下のような処理を任意の言語（python,R）で実装すればいいでしょう。

1.元の説明変数と被説明変数の関係を行列式を用いてY=Xに対して線対称の値に変換する。
2.Logit関数を前提としてモデルパラメータを計算する。
3.理論値を取得し、改めてY=Xに対して線対称の値に変換する。

pythonの場合は以下のような感じになります。

def func1(x) :   #### toy data 作成用関数
    y = 3.89*x**3 - 5.56*x**2 + 2.65*x +0.013 + stats.norm.rvs(loc=0,scale=0.05,size=len(x))
    return y

x = np.linspace(0,1,11)
y = func1(x)

def reflection(x, m):    #### y=mxに対して線対称の値に変換する関数
    a = np.array([[1 - m**2, 2*m],
                  [2*m, m**2 - 1]])
    a = a / (1 + m**2)
    ax = a@x
    return ax

x_translate = reflection(np.array([x,y]),1)

#### Logitモデルのパラメータ推計（リンク関数はY=aX+bと仮定）

import statsmodels.api as sm
X = sm.add_constant(x_translate[0,:])
logit_model = sm.Logit(x_translate[1,:],X)
result = logit_model.fit()


#### 結果を再び変換
res = reflection(np.array([x_translate[0,:],result.predict()]),1)

Comments

[yumetodo]

Logit関数について自分なりに調べたあとでBA決定します。

Answer 2

y = ax^3+bx^2+(1-a-b)x=a(x^3-x) + b(x^2-x)+ x として最小二乗法。

Comments

[yumetodo]

そもそも式変形すれば一本の式にできることに頭が回ってなかったです。
実際に試してみてからBA決定します。