2017-09-11 1 views
0

J'essaie d'utiliser constrOptim pour optimiser la somme des erreurs carrées d'une régression multiple linéaire. L'équation principale devrait être D = Beta1*Xa+Beta2*Xb+Beta3*Xc+Beta4*Xd, avec D,Xa,Xb,Xc,Xd à partir d'un fichier .csv importé, et les Beta s sont les coefficients que je veux trouver, en minimisant les erreurs quadratiques.Mauvais résultat de la fonction constroptim

Jusqu'à présent, j'importé le file.csv à R, nommé chaque colonne comme Ds, Xa, Xb, Xc, Xd, a créé le objfunction = function(Beta1,Beta2,Beta3,Beta4)'sum(E²)'=(sum(D) - sum(Beta1*Xa+Beta2*Xb+Beta3*Xc+Beta4*Xd))^2)

créé la matrice 'C' et le vecteur 'd' pour configurer les contraintes qui devraient limiter les bêta à < = 0. Je ne sais pas comment trouver la région possible, bien que j'ai utilisé les valeurs initiales qui ont fait fonctionner la fonction.

Voici le code:

> Tabela= read.table("Simulacao.csv", header=T, sep= ";") 
> Tabela 
    D A B C D.1 
1 -1 1 -1 0 0 
2 4 0 0 1 -1 
3 4 1 0 -1 0 
4 0 0 1 0 -1 
5 -2 1 0 0 -1 
> Ds= Tabela[,1] 
> Xa= Tabela[,2] 
> Xb= Tabela[,3] 
> Xc= Tabela[,4] 
> Xd= Tabela[,5] 
> simulaf= function(x1,x2,x3,x4) { 
+ Ds= Tabela[,1] 
+ Xa= Tabela[,2] 
+ Xb= Tabela[,3] 
+ Xc= Tabela[,4] 
+ Xd= Tabela[,5] 
+ J=sum(Ds) 
+ H=sum(x1*Xa+x2*Xb+x3*Xc+x4*Xd) 
+ sx=(J-H)^2 
+ return(sx) 
+ } 
> s= function(x) {simulaf(x[1],x[2],x[3],x[4])} 
> d= c(0,0,0,0) 
> C= matrix(c(-1,0,0,0,0,-1,0,0,0,0,-1,0,0,0,0,-1),nrow=4,ncol=4,byrow=T) 
> constrOptim(c(-1,-1,-1,-1),s,NULL,C,d) 
$par 
[1] -0.2608199 -0.8981110 -1.1095961 -1.9274866 

Le résultat devrait être que je pense:

$par 
[1] -0.125 0 -0.5 -0.875 

Après des recherches, mes conclusions sont que cela pourrait être parce que je suis en utilisant de mauvaises valeurs initiales, problème de paramétrage (ne comprends pas pourquoi c'est nécessaire) ou si c'est simplement que je l'ai mal programmé.

Que dois-je faire pour résoudre ce problème?

Répondre

1

La formule de la somme des carrés des écarts est

sum((y - yhat)^2) 

et pas

(sum(y) - sum(yhat))^2 

yhat est la valeur prédite. En outre, si vos seules contraintes sont que les bêtas estimés devraient être négatifs (ce qui est un peu bizarre, généralement vous voulez qu'ils soient positifs mais peu importe), alors vous n'avez pas besoin de constrOptim. La norme ou nlminb fonctionnera avec des contraintes de boîte.