OLSのロバスト標準誤差

import numpy as np
import matplotlib.pyplot as plt

n = 100
np.random.seed(0)
x1 = np.random.uniform(0, 5, size=(n, 1))
x0 = np.ones(shape=(n, 1))
X = np.append(x0, x1, axis=1)
beta = np.array([3, 5])

# 均一分散の場合
e = np.random.normal(loc=0, scale=1, size=n)
y = X @ beta + e

# plot
fig, ax = plt.subplots()
ax.scatter(x1, y)

<matplotlib.collections.PathCollection at 0x7f6b38343a00>

beta_hat = np.linalg.inv(X.T @ X) @ X.T @ y
beta_hat

array([3.22215108, 4.987387  ])

X.T @ X

array([[100.        , 236.39691976],
       [236.39691976, 766.62957534]])

(X[0, ] @ X[0,].T)**(-1)

0.11723457858134044