Number of Parameters, AIC, and BIC
Simple Simulations
Brunton and Kutz (2019) explained that error increases as number of parameters increases beyond some value. But, in this case, error is measured from the true model, which is unobservable. As is shown at this website, error measured from the observed data decreases as number of parameters increases. Hence, it is difficult to select a model as true one based on information of errors. Brunton and Kutz (2019) explained usage of information criteria, AIC and BIC, in selecting the mode. I wrote simple simulation scripts to study on AIC and BIC, in the framework of Code 4.14 of Brunton and Kutz (2019).
True model:
Data
Models used
Polynomial model of degree p
Information criteria (Konishi and Kitagawa, 2008)
Two simulation scripts were prepared for the above true model.
(1) For fixed data, models f_pred() were compared.
(2) For sample sizes, n_pnts=20, 100, 1000, information criteria were compared.
The scripts were run in Python 3.12.4/Anaconda. Programs downloaded from https://www.python.org/ caused error(2024.06.26).
Script for the case of (1) is shown in Listing 1.
Run the script in Listing 1, graphs in Figure 1 are shown.
Figure 1
Information criteria, AIC and BIC, show trends that is minimum at the value of the true model.
Close the window of Figure 1, the next graph is presented (Figure 2).
Figure 2
The histogram shows frequencies of selections of models with respective polynomial degrees. Selections by information criteria sometimes lead into errors.
The script to study effects by sample sizes is shown in Listing 2.
Run the script, then the graph in Figure 3 is displayed.
Figure 3
Frequencies of selected models are shown by line graphs to draw the results of the six conditions clearly. With larger sample size, frequency of correct selection increases, although differences between sample sizes 100 and 1000 are negligible.
import math
import numpy as np
import matplotlib.pyplot as plt
M = 5
n = 1000
L = 4
n_pnts = 100
x = np.linspace(0, L, n_pnts)
f = x**2
# True value
AICs = np.empty((n, M+1))
BICs = np.empty((n, M+1))
aic_choice_p = []
bic_choice_p = []
for i_trial in range(n):
# True value + Noise
fn = f + 0.1 *
np.random.normal(size=n_pnts)
for p in range(M+1):
A =
np.empty((n_pnts, p+1))
for
j in range(p+1):
A[:,j] = x ** j
b =
np.linalg.pinv(A) @ fn
f_pred = A @ b
sse
= np.sum((fn - f_pred)**2)
sgm
= (sse / n_pnts) ** 0.5
max_log_L = -0.5 * n_pnts * math.log(2*math.pi) \
- n_pnts * math.log(sgm) \
-sse /(2 * (sgm**2))
AICs[i_trial][p] = -2 * max_log_L + 2 * (p + 1 + 1)
BICs[i_trial][p] = -2 * max_log_L + (p +
1 + 1) * math.log(n_pnts)
c_aic =
np.argmin(AICs[i_trial])
aic_choice_p.append(c_aic)
c_bic =
np.argmin(BICs[i_trial])
bic_choice_p.append(c_bic)
plt.figure(figsize=(12,5))
plt.subplot(1,2,1)
Mean_AICs = np.mean(AICs, axis=0)
plt.plot(range(M+1), Mean_AICs, '^',
c='r', lw=0, label='AIC')
Mean_BICs = np.mean(BICs, axis=0)
plt.plot(range(M+1), Mean_BICs, 'v', c
= 'g', lw=0, label='BIC')
plt.boxplot(AICs,
positions=np.arange(0,M+1))
plt.boxplot(BICs,
positions=np.arange(0,M+1))
plt.xticks(np.arange(0, M+1))
plt.xlabel('Polynomial
Degree',fontsize=14)
plt.ylabel('Information Criterion',
fontsize=14)
plt.title('Degree >= 0',
fontsize=16)
plt.legend(fontsize=14)
plt.subplot(1,2,2)
plt.plot(np.arange(M+1)[2:],
Mean_AICs[2:], '^-', c='r', lw=1, label='AIC')
plt.plot(np.arange(M+1)[2:],
Mean_BICs[2:], 'v-', c = 'g', lw=1, label='BIC')
plt.xticks(np.arange(2, M+1))
plt.xlabel('Polynomial Degree',
fontsize=14)
plt.ylabel('Information Criterion',
fontsize=14)
plt.boxplot(AICs[:,2:],
positions=np.arange(2, M+1), boxprops={'color':'r'},
flierprops={'markerfacecolor':'r'})
plt.boxplot(BICs[:,2:],
positions=np.arange(2, M+1), boxprops={'color':'g'},
flierprops={'markerfacecolor':'g'})
plt.title('Degree >= 2',
fontsize=16)
plt.legend()
plt.tight_layout()
plt.show()
plt.title('Model Selection',
fontsize=16)
plt.hist([aic_choice_p, bic_choice_p],
bins=np.arange(-0.5, M+1, 1.0),
color=['r', 'g'], label=['AIC', 'BIC'])
plt.ylabel('Frequency', fontsize=14)
plt.xlabel('Degree', fontsize=14)
plt.legend(fontsize=14)
plt.tight_layout()
plt.show()
Listing 2
import math
import numpy as np
import matplotlib.pyplot as plt
M = 10
n = 1000
L = 4
choice_aic = np.zeros((3, M+1),
dtype=int)
choice_bic = np.zeros((3, M+1),
dtype=int)
n_pnts_L = [20, 100, 1000]
for i, n_pnts in enumerate(n_pnts_L):
x = np.linspace(0, L,
n_pnts)
f = x**2
# True value
AICs = np.empty((n, M+1))
BICs = np.empty((n, M+1))
for i_trial in range(n):
# True value +
Noise
fn =
f + 0.1 * np.random.normal(size=n_pnts)
for
p in range(M+1):
A = np.empty((n_pnts, p+1))
for j in range(p+1):
A[:,j] = x ** j
b
= np.linalg.pinv(A) @ fn
f_pred = A @ b
sse = np.sum((fn - f_pred)**2)
sgm = (sse / n_pnts) ** 0.5
max_log_L = -0.5 * n_pnts * math.log(2*math.pi) \
- n_pnts * math.log(sgm) \
-sse /(2 * (sgm**2))
AICs[i_trial][p] = -2 * max_log_L + 2 * (p + 1 + 1)
BICs[i_trial][p] = -2 * max_log_L + (p + 1 + 1) * math.log(n_pnts)
c_aic = np.argmin(AICs[i_trial])
choice_aic[i, c_aic] += 1
c_bic = np.argmin(BICs[i_trial])
choice_bic[i, c_bic] += 1
for i in range(3):
print()
print(f'sample size =
{n_pnts_L[i]}')
print(f'{"Degree":
<15s}', end='')
for j in range(M+1):
print(f'{j:>10d}', end='')
print()
print(f'{"Freq.(AIC)":<15s}', end='')
for j in range(M+1):
print(f'{choice_aic[i,j]:>10d}', end='')
print()
print(f'{"Freq.(BIC)":<15s}', end='')
for j in range(M+1):
print(f'{choice_bic[i,j]:>10d}', end='')
print()
lines = ['--', ':', '-']
for i in range(3):
plt.plot(range(M+1),
choice_aic[i], ls=lines[i], c='r',
label=f'AIC-{n_pnts_L[i]}')
plt.plot(range(M+1),
choice_bic[i], ls=lines[i], c='g',
label=f'BIC-{n_pnts_L[i]}')
plt.title('Model Selection',
fontsize=16)
plt.ylabel('Frequency', fontsize=14)
plt.xlabel('Degree', fontsize=14)
plt.legend()
plt.tight_layout()
plt.show()
Brunton, S. L., and Kutz, J. N. (2019). Data-driven science and engineering: Machine learning, dynamic systems, and control. Cambridge University Press.
Konishi, S., and Kitagawa, G. (2008). Information criteria and statistical modeling. Springer.