{"word_count": 451, "character_count": 4375}
Pour travailler un peu la manipulation de matrices, dans le cadre du cours d'analyse de données, quelques lignes de code.La matrice de données est ici&#62; BD=read.table("https://perso.univ-rennes1.fr/arthur.charpentier/securite.txt",header=TRUE)&#62; M=as.matrix(BD[,2:10])&#62; colnames(M)=names(BD)[2:10]&#62; rownames(M)=BD$dep&#62; head(M[,1:5])    infra  vols  eco crim  vmaD1  44.11 29.51 4.14 2.79 0.08D2  45.97 26.29 3.97 4.54 0.07D3  38.83 21.90 4.61 3.33 0.03D4  49.68 30.11 5.28 4.54 0.07D5  47.67 27.19 4.94 3.91 0.05D6 109.21 77.71 6.87 5.25 0.24La décomposition spectrale de&#62; round(t(M)%*%M)       infra   vols   eco  crim vma  vvi  camb  infra 308599 191433 29179 21290 642 7042 30951 vols  191433 119961 17859 12967 404 4427 19335 eco    29179  17859  2946  2042  59  648  2942  crim   21290  12967  2042  1564  42  448  2129  vma      642    404    59    42   2   19    62   vvi     7042   4427   648   448  19  237   665  camb   30951  19335  2942  2129  62  665  3200  roul   55363  34745  5089  3792 113 1217  5612 auto   24841  15747  2214  1676  53  545  2520(j'ai viré les deux dernières colonnes pour des raisons de présentation) donne ici&#62; eigen(t(M)%*%M)$values[1] 4.471766e+05 1.239465e+03 3.462936e+02 1.764038e+02 1.303316e+02 6.623592e+01 2.624739e+01 1.761801e+01 1.564292e-01 $vectors(etc). Pour rappel, les liens entre la diagonalisation et l'analyse spectrale se voit avec l'écriture&#62; P=eigen(t(M)%*%M)$vectors&#62; P%*%diag(eigen(t(M)%*%M)$values)%*%t(P)Maintenant, deux remarques. Tout d'abord le lien entre la décomposition spectrale de cette matrice, et la décomposition en valeurs singulières&#62; sqrt(eigen(t(M)%*%M)$values)[1] 668.7126681  35.2060415  18.6089656  13.2817097  11.4162872   8.1385454   5.1232210   4.1973809   0.3955113&#62; svd(M)$d[1] 668.7126681  35.2060415  18.6089656  13.2817097  11.4162872   8.1385454   5.1232210   4.1973809   0.3955113 $u             [,1]          [,2]          [,3]          [,4]         [,5]          [,6]         [,7]         [,8]          [,9] [1,] -0.08156446 -0.0632514738 -2.596265e-02  0.0960587227 -0.064759395 -0.0430947134 -0.024719699 -0.074347281  0.0550356487 [2,] -0.08087892  0.0541859480 -5.777492e-02 -0.0628249751 -0.099084701 -0.1311683858  0.083484434  0.011872601 -0.0195034904(etc) et la décomspotion spectrale de l'autre produit matriciel&#62; sqrt(eigen(M%*%t(M))$values) [1] 6.687127e+02 3.520604e+01 1.860897e+01 1.328171e+01 1.141629e+01 8.138545e+00 5.123221e+00 4.197381e+00 3.955113e-01 7.546505e-06 6.656212e-06 5.285131e-06 4.270035e-06(etc).Maintenant, pour comprendre un peu mieux la recherche de composantes pricinpales, considérons les deux variables&#62; sM=M[,c(1,3)]&#62; plot(sM)&#160;On a tout intérêt à centrer et réduire les variables (ou changer la métrique, comme vu en cours)&#62; sMcr=sM&#62; for(j in 1:2) sMcr[,j]=(sMcr[,j]-mean(sMcr[,j]))/sd(sMcr[,j])&#62; plot(sMcr)Avant de faire des projections sur un axe, introduisons deux fonctions&#62; prod_scal=function(x,u){+   ps=rep(NA,nrow(x))+   for(i in 1:nrow(x)) ps[i]=sum(x[i,]*u)+   return(ps)+ }&#62; proj=function(x,u){+   px=x+   for(j in 1:length(u)){+     px[,j]=prod_scal(x,u)/sqrt(sum(u^2))*u[j]  +   }+   return(px)+ }Par exemple, si on projette sur l'axe des abscisses,&#62; points(proj(sMcr,c(1,0)),col="blue")On peut ensuite chercher la direction de l'axe qui nous donnera le nuage de points qui aura la plus grande inertie&#62; inertie=function(x) sum(x^2)&#62; Theta=seq(0,3.141592657,length=101)&#62; V=unlist(lapply(Theta,function(theta) inertie(proj(sMcr,c(cos(theta),sin(theta))))))&#62; plot(Theta,V,type='l')&#62; (angle=optim(0,function(theta) -inertie(proj(sMcr,c(cos(theta),sin(theta)))))$par)[1] 0.7853516Visuellement, on obtient&#62; plot(sMcr)&#62; psMcr=proj(sMcr,c(cos(angle),sin(angle)))&#62; points(psMcr,col="blue")On notera que l'axe qui donne l'inertie maximale est lié aux vecteurs propres de la décomposition spectrale (celui associé à la plus grande valeur propre)&#62; c(cos(angle),sin(angle))[1] 0.7071397 0.7070738&#62; eigen(t(sMcr)%*%sMcr)$values[1] 159.64663  26.35337 $vectors          [,1]       [,2][1,] 0.7071068 -0.7071068[2,] 0.7071068  0.7071068On continuera demain, à manipuler des matrices de données, de faire des projections, avant de commencer l'analyse en composantes principales, à proprement parler, la semaine prochaine