Pourquoi existe-t-il de grands coefficients pour les polynômes d'ordre supérieur

13

Dans le livre de Bishop sur l'apprentissage automatique, il traite du problème de l'ajustement d'une courbe d'une fonction polynomiale à un ensemble de points de données.

Soit M l'ordre du polynôme ajusté. Il déclare que

Nous voyons qu'à mesure que M augmente, l'amplitude des coefficients augmente généralement. En particulier pour le polynôme M = 9, les coefficients sont devenus finement ajustés aux données en développant de grandes valeurs positives et négatives de sorte que la fonction polynomiale correspondante correspond exactement à chacun des points de données, mais entre les points de données (en particulier près des extrémités des gamme) la fonction présente les grandes oscillations.

Je ne comprends pas pourquoi de grandes valeurs impliquent de mieux ajuster les points de données. Je pense que les valeurs deviendraient plus précises après la décimale pour un meilleur ajustement.

Abhishek Bhatia
la source
le livre considère x à 10 points également espacés dans [ 0 , 1 ]ϵ est gaussien avec une moyenne nulle et une «petite» variance (en considérant donc l'ajustement d'un polynôme à 9 dimensions avec 10 points ...y=sjen(2πX)+ϵ[0,1]ϵ
seanv507

Réponses:

18

Il s'agit d'un problème bien connu avec les polynômes d'ordre élevé, connu sous le nom de phénomène de Runge . Numériquement, il est associé à un mauvais conditionnement de la matrice de Vandermonde , ce qui rend les coefficients très sensibles aux petites variations des données et / ou à l'arrondi des calculs (ie le modèle n'est pas identifiable de manière stable ). Voir aussi cette réponse sur le SciComp SE.

Il existe de nombreuses solutions à ce problème, par exemple l' approximation de Tchebychev , le lissage des splines et la régularisation de Tikhonov . La régularisation de Tikhonov est une généralisation de la régression des crêtes , pénalisant une norme du vecteur de coefficient θ , où pour lisser la matrice de poids Λ est un opérateur dérivé. Pour pénaliser les oscillations, vous pouvez utiliser Λ θ = p [ x ] , où p [ x ]||Λθ]||θΛΛθ=p[X]p[X] est le polynôme évalué au niveau des données.

EDIT: La réponse de l'utilisateur hxd1011 note que certains des problèmes de mauvais conditionnement numérique peuvent être résolus en utilisant des polynômes orthogonaux, ce qui est un bon point. Je voudrais cependant noter que les problèmes d'identification des polynômes d'ordre élevé demeurent. Autrement dit, le mauvais conditionnement numérique est associé à la sensibilité aux perturbations "infinitésimales" (par exemple arrondi), tandis que le mauvais conditionnement "statistique" concerne la sensibilité aux perturbations "finies" (par exemple les valeurs aberrantes; le problème inverse est mal posé ).

Les méthodes mentionnées dans mon deuxième paragraphe concernent cette sensibilité aberrante . Vous pouvez considérer cette sensibilité comme une violation du modèle de régression linéaire standard, qui en utilisant une inadéquation suppose implicitement que les données sont gaussiennes. Splines et régularisation de Tikhonov gèrent cette sensibilité aberrante en imposant une finesse préalable à l'ajustement. L'approximation de Chebyshev traite cela en utilisant un inadéquat L appliqué sur le domaine continu , c'est-à-dire pas seulement aux points de données. Bien que les polynômes de Chebyshev soient orthogonaux (par rapport à un certain produit intérieur pondéré), je crois que s'ils étaient utilisés avec un L 2 inadapté aux données, ils seraient toujoursL2LL2 ont une sensibilité aux valeurs aberrantes.

GeoMatt22
la source
@ hxd1011 c'est vrai, et j'ai donné l'exemple des polynômes de Chebyshev. Est-ce que les polynômes orthogonaux résoudront toujours le problème dans la pratique, s'il y a des valeurs aberrantes et que l'inadéquation des données est toujours ? Je pense que le phénomène Runge causerait toujours des problèmes de fiabilité dans les coefficients dans ce cas (c'est-à-dire pour les perturbations finies / grandes sur les données, par rapport à l'arrondi infinitésimal.)L2
GeoMatt22
1
p
1
Où puis-je en savoir plus sur ce métier de "mauvais conditionnement de la matrice vandermonde"?
Matthew Drury
@MatthewDrury Je fais habituellement aussi l'approche empirique suggérée par hxd1011. Cependant, après votre requête, un rapide Google a révélé un article récent qui pourrait également être intéressant: à quel point les matrices Vandermonde sont-elles mauvaises? (VY Pan, 2015) . (Par exemple, il explique pourquoi les matrices DFT sont Vandermonde mais pas mal conditionnées.)
GeoMatt22
Merci @ GeoMatt22. Désolé de vous avoir fait google pour moi, ai-je demandé car je pensais que vous pouviez avoir des sources personnelles préférées.
Matthew Drury
8

La première chose que vous voulez vérifier, c'est si l'auteur parle de polynômes bruts par rapport à des polynômes orthogonaux .

Pour les polynômes orthogonaux. le coefficient ne devient pas "plus grand".

Voici deux exemples d'expansion polynomiale d'ordre 2 et 15. D'abord, nous montrons le coefficient d'expansion du second ordre.

summary(lm(mpg~poly(wt,2),mtcars))

Call:
lm(formula = mpg ~ poly(wt, 2), data = mtcars)

Residuals:
   Min     1Q Median     3Q    Max 
-3.483 -1.998 -0.773  1.462  6.238 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)   20.0906     0.4686  42.877  < 2e-16 ***
poly(wt, 2)1 -29.1157     2.6506 -10.985 7.52e-12 ***
poly(wt, 2)2   8.6358     2.6506   3.258  0.00286 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.651 on 29 degrees of freedom
Multiple R-squared:  0.8191,    Adjusted R-squared:  0.8066 
F-statistic: 65.64 on 2 and 29 DF,  p-value: 1.715e-11

Ensuite, nous montrons le 15ème ordre.

summary(lm(mpg~poly(wt,15),mtcars))

Call:
lm(formula = mpg ~ poly(wt, 15), data = mtcars)

Residuals:
    Min      1Q  Median      3Q     Max 
-5.3233 -0.4641  0.0072  0.6401  4.0394 

Coefficients:
               Estimate Std. Error t value Pr(>|t|)    
(Intercept)     20.0906     0.4551  44.147  < 2e-16 ***
poly(wt, 15)1  -29.1157     2.5743 -11.310 4.83e-09 ***
poly(wt, 15)2    8.6358     2.5743   3.355  0.00403 ** 
poly(wt, 15)3    0.2749     2.5743   0.107  0.91629    
poly(wt, 15)4   -1.7891     2.5743  -0.695  0.49705    
poly(wt, 15)5    1.8797     2.5743   0.730  0.47584    
poly(wt, 15)6   -2.8354     2.5743  -1.101  0.28702    
poly(wt, 15)7    2.5613     2.5743   0.995  0.33459    
poly(wt, 15)8    1.5772     2.5743   0.613  0.54872    
poly(wt, 15)9   -5.2412     2.5743  -2.036  0.05866 .  
poly(wt, 15)10  -2.4959     2.5743  -0.970  0.34672    
poly(wt, 15)11   2.5007     2.5743   0.971  0.34580    
poly(wt, 15)12   2.4263     2.5743   0.942  0.35996    
poly(wt, 15)13  -2.0134     2.5743  -0.782  0.44559    
poly(wt, 15)14   3.3994     2.5743   1.320  0.20525    
poly(wt, 15)15  -3.5161     2.5743  -1.366  0.19089    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.574 on 16 degrees of freedom
Multiple R-squared:  0.9058,    Adjusted R-squared:  0.8176 
F-statistic: 10.26 on 15 and 16 DF,  p-value: 1.558e-05

Notez que nous utilisons des polynômes orthogonaux , donc le coefficient de l'ordre inférieur est exactement le même que les termes correspondants dans les résultats de l'ordre supérieur. Par exemple, l'ordonnée à l'origine et le coefficient pour le premier ordre sont 20,09 et -29,11 pour les deux modèles.

dix6

> summary(lm(mpg~poly(wt,15, raw=T),mtcars))

Call:
lm(formula = mpg ~ poly(wt, 15, raw = T), data = mtcars)

Residuals:
    Min      1Q  Median      3Q     Max 
-5.6217 -0.7544  0.0306  1.1678  5.4308 

Coefficients: (3 not defined because of singularities)
                          Estimate Std. Error t value Pr(>|t|)
(Intercept)              6.287e+05  9.991e+05   0.629    0.537
poly(wt, 15, raw = T)1  -2.713e+06  4.195e+06  -0.647    0.526
poly(wt, 15, raw = T)2   5.246e+06  7.893e+06   0.665    0.514
poly(wt, 15, raw = T)3  -6.001e+06  8.784e+06  -0.683    0.503
poly(wt, 15, raw = T)4   4.512e+06  6.427e+06   0.702    0.491
poly(wt, 15, raw = T)5  -2.340e+06  3.246e+06  -0.721    0.480
poly(wt, 15, raw = T)6   8.537e+05  1.154e+06   0.740    0.468
poly(wt, 15, raw = T)7  -2.184e+05  2.880e+05  -0.758    0.458
poly(wt, 15, raw = T)8   3.809e+04  4.910e+04   0.776    0.447
poly(wt, 15, raw = T)9  -4.212e+03  5.314e+03  -0.793    0.438
poly(wt, 15, raw = T)10  2.382e+02  2.947e+02   0.809    0.429
poly(wt, 15, raw = T)11         NA         NA      NA       NA
poly(wt, 15, raw = T)12 -5.642e-01  6.742e-01  -0.837    0.413
poly(wt, 15, raw = T)13         NA         NA      NA       NA
poly(wt, 15, raw = T)14         NA         NA      NA       NA
poly(wt, 15, raw = T)15  1.259e-04  1.447e-04   0.870    0.395

Residual standard error: 2.659 on 19 degrees of freedom
Multiple R-squared:  0.8807,    Adjusted R-squared:  0.8053 
F-statistic: 11.68 on 12 and 19 DF,  p-value: 2.362e-06
Haitao Du
la source
Je ne suis pas sûr que la syntaxe soit correcte, mais pourquoi ne contrastez-vous pas les résultats de v orthogonal brut avec quelque chose commesummary(lm(mpg~poly(wt,2),mtcars)); summary(lm(mpg~poly(wt,5),mtcars)); summary(lm(mpg~ wt + I(wt^2),mtcars)); summary(lm(mpg~ wt + I(wt^2) + I(wt^3) + I(wt^4) + I(wt^5),mtcars))
Antoni Parellada
@AntoniParellada bonne suggestion, je vais réviser. BTW, nous pouvons facilement faire une expansion brute parpoly(x,2,raw=T)
Haitao Du
Belle astuce ... Je suppose que vous pouvez vous en tenir à 15 et le faire summary(lm(mpg~poly(wt,15, raw=T),mtcars)). Effet massif dans les coefficients!
Antoni Parellada
Un commentaire sur ma réponse de @ seanv507 m'a rendu curieux de savoir ce qui suit. Si vous utilisez des polynômes orthogonaux et souhaitez réduire la sensibilité aux valeurs aberrantes, une régression de crête standard serait-elle alors suffisante? Ou les polynômes d'ordre supérieur, plus oscillatoires, nécessiteraient-ils encore des poids ~ ordre? (Je pense que ce dernier, comme par exemple une matrice DFT est orthogonale, mais la sous-pondération des hautes fréquences serait toujours nécessaire. J'ai eu une expérience (désagréable) avec ce cas particulier!)
GeoMatt22
3

Abhishek, vous avez raison, l'amélioration de la précision des coefficients améliorera la précision.

Nous voyons qu'à mesure que M augmente, l'amplitude des coefficients augmente généralement. En particulier pour le polynôme M = 9, les coefficients sont devenus finement ajustés aux données en développant de grandes valeurs positives et négatives de sorte que la fonction polynomiale correspondante correspond exactement à chacun des points de données, mais entre les points de données (en particulier près des extrémités des plage) la fonction présente de grandes oscillations.

Je pense que le problème de l'ampleur n'est pas pertinent pour le point de vue de Bishop - que l'utilisation d'un modèle compliqué sur des données limitées conduit à un `` sur-ajustement ''. Dans son exemple, 10 points de données sont utilisés pour estimer un polynôme à 9 dimensions (c'est-à-dire 10 variables et 10 inconnues).

Si nous ajustons une onde sinusoïdale (pas de bruit), alors l'ajustement fonctionne parfaitement, car les ondes sinusoïdales [sur un intervalle fixe] peuvent être approximées avec une précision arbitraire en utilisant des polynômes. Cependant, dans l'exemple de Bishop, nous avons une certaine quantité de «bruit» que nous ne devrions pas adapter. Pour ce faire, nous maintenons le nombre de points de données par rapport au nombre de variables de modèle (coefficients polynomiaux) ou en utilisant la régularisation. Ajustement polynomial du 9ème ordre à 10 points de datation sur (0,1)

La régularisation impose des contraintes «douces» au modèle (par exemple, dans la régression de crête), la fonction de coût que vous essayez de minimiser est une combinaison d '«erreur d'ajustement» et de complexité du modèle: par exemple, dans la régression de crête, la complexité est mesurée par la somme des coefficients carrés. cela entraîne un coût sur la réduction de l'erreur - l'augmentation des coefficients ne sera autorisée que si elle a une réduction suffisamment grande de l'erreur d'ajustement [la taille de cette erreur est suffisamment grande est spécifiée par un multiplicateur sur le terme de complexité du modèle]. Par conséquent, l'espoir est qu'en choisissant le multiplicateur approprié, nous ne nous adapterons pas à un petit terme de bruit supplémentaire, car l'amélioration de l'ajustement ne justifie pas l'augmentation des coefficients.

Vous avez demandé pourquoi des coefficients élevés améliorent la qualité de l'ajustement. Essentiellement, la raison en est que la fonction estimée (sin + bruit) n'est pas un polynôme, et les grands changements de courbure nécessaires pour approximer l'effet de bruit avec des polynômes nécessitent de grands coefficients.

Notez que l'utilisation de polynômes orthogonaux n'a aucun effet (j'ai ajouté un décalage de 0,1 juste pour que les polynômes orthogonaux et bruts ne soient pas superposés)

require (penalized)
poly_order<-9
x_long<-seq(0,1, length.out = 100)
nx<-10
x<-seq(0,1, length.out = nx)
noise<- rnorm(nx, 0, 1)
noise_scale<-0.2
y<-sin(2*pi*x)+noise_scale*noise

training_data<-data.frame(x=x,y=y)
y_long<-sin(2*pi*x_long)

plot(x,y, col ='blue',ylim=c(-1.5,1.5))
lines(x_long,y_long,col='green')

polyfit_raw<-lm(y~poly(x,poly_order,raw=TRUE),data=training_data)
summary(polyfit_raw)

polyfit_raw_ridge1<-penalized(y,~poly(x,poly_order,raw=TRUE), model='linear', data=training_data, lambda2=0.0001, maxiter=10000, standardize=TRUE)

polyfit_orthog<-lm(y~poly(x,poly_order),data=training_data)
summary(polyfit_orthog)

pred_raw<-predict(polyfit_raw,data.frame(x=x_long))
pred_ortho<-predict(polyfit_orthog,data.frame(x=x_long))
pred_raw_ridge<-predict(polyfit_raw_ridge1,data=data.frame(x=x_long))[,'mu']
lines(x_long,pred_raw,col='red')
# add 0.1 offset to make visible
lines(x_long,pred_ortho+0.1,col='black')
lines(x_long,pred_raw_ridge,col='purple')
legend("bottomleft",legend=c('data sin(2 pi x) + noise','sin(2 pi x)', 
                             'raw poly','orthog poly +0.1 offset','raw poly + ridge regression'),
       fill=c('blue','green','red','black','purple'))
seanv507
la source