Monétisation et marchandisation des données personnelles

Pour le dernier billet de la série autour du big data sur variances.eu, Fabrice Rochelandet revient sur la monétisation et la marchandisation des données personnelles,

L’exploitation marchande des données personnelles (DP) a donné lieu à une importante littérature en économie et en marketing (Hui & Png [2006] ; Anderson & Moore [2006] ; Rochelandet [2010], Tucker [2012], Acquisti et al. [2016]). La collecte et la transformation de ces ressources informationnelles par les entreprises y sont présentées comme une source de bien-être pour les consommateurs, d’innovations et d’opportunités commerciales inédites, mais également de menaces sur la vie privée et d’externalités négatives pécuniaires. Dans ce dernier cas, les exploitants de DP peuvent par exemple opérer des stratégies de discrimination par les prix sur les marchés de biens de consommation courants, de l’assurance, des prêts bancaires, du travail… en fonction du profilage des individus : les individus les plus riches, les moins solvables ou présentant des risques de santé paieront un prix supérieur.

Les entreprises collectent des données sur leurs clients pour les exploiter selon des modalités très différentes : les identifier, communiquer ou interagir avec eux de façon plus ou moins ciblée, différencier leurs prix ou la qualité de leurs offres en fonction de leurs profils (marché du crédit bancaire par exemple), animer une communauté en ligne, faciliter l’appariement entre les utilisateurs de leurs services (services de réseautage social ou de rencontre en ligne par exemple), vendre ces données. Les moyens de collecte et de production de ces informations sont eux-mêmes très diversifiés. Les individus peuvent ainsi divulguer volontairement leurs données par déclaration (formulaires d’abonnement et d’inscription, données de transaction, etc.) ou en acceptant que leurs données soient enregistrées (traces de navigation, localisation géographique, mots clés extraits de leurs commentaires et discussions en ligne, etc.). L’exposition de soi ex post (Rallet & Rochelandet [2011]) s’avère également une source majeure d’informations sur l’internet via la collecte massive de données enregistrées par des traceurs et leur traitement par des algorithmes permettant de déduire des informations à valeur ajoutée à partir de données brutes très disparates.

Dans certains cas, la collecte des données d’usage bénéficie essentiellement à leurs exploitants. Par exemple, les sites de vente en ligne analysent l’audience web et son impact sur les ventes. La connaissance des pages consultées antérieurement par les internautes leur révèle des informations précieuses pour optimiser leur campagne promotionnelle. De même, les avis postés en ligne permettent d’alimenter au moindre coût les systèmes de recommandation de plateformes comme Amazon et Netflix qui orientent les choix de leurs clients, les fidélisent et différencient leurs services.

[à suivre…]

Acquisti, A., C.R. Taylor, L. Wagman (2016) “The economics of privacy”, Journal of Economic Literature, vol. 54, n°2, pp.442-492.
Anderson R., T. Moore (2006) “The economics of security”, Science, vol.314, n°5799, pp.610–613.
Brandimarte, L., A. Acquisti, G. Loewenstein (2011) “Misplaced confidences: Privacy and the control paradox”, Social Psychological and Personality Science,
Cecere, G. F. Le Guel, F. Rochelandet (2015) “Les modèles d’affaires numériques sont-ils trop indiscrets ? Une analyse empirique”, Réseaux, vol.33, n°189, pp.77-101.
Hui, K.L., I.P.L. Png (2006) “The economics of privacy”, in Economics and Information Systems, Handbooks in Information Systems, vol.1, chapitre 9, ed. Terrence Hendershott, Elsevier.
Pucheral, P., A. Rallet, F. Rochelandet, C. Zolinski (2016) “La Privacy by design : une fausse bonne solution aux problèmes de protection des données personnelles soulevés par l’Open data et les objets connectés ?”, Légicom, 2016/1, n°56, pp.89-99.
Rallet, A., F. Rochelandet (2011) “La régulation des données personnelles face au web relationnel : une voie sans issue ?”, Réseaux, vol.29, n°167, pp.19-47.
Rochelandet, F. (2010) Economie des données personnelles et de la vie privée, La Découverte.
Goldfarb, A., C. Tucker (2011) “Online display advertising: Targeting and obtrusiveness”, Marketing Science, vol.30, pp.389-404.
Tucker, C. E. (2012) “The economics of advertising and privacy”, International journal of Industrial organization, vol.30, n°3, pp.326–329.
White, T., D. Zahay, H. Thorbjornsen, S. Shavitt (2008) “Getting too personal: Reactance to highly personalized email solicitations”, Marketing Letters, vol.19, n°1, pp.39–50.

Assurance, solidarité ou individualisation

La semaine dernière, à l’ensae, on a eu une introduction générale à la tarification dans une population hétérogène et la segmentation. Dans deux semaines, je participerais à un colloque à la Faculté de Droit d’Amiens sur l’utilisation des données personnelles en assurance. Je me disais que ça pourrait être l’occasion de (re)mettre par écrit certaines réflexions générales sur le sujet…

La segmentation en matière d’assurance évoque la classification qu’un assureur opère pour fixer le niveau des primes, de telle sorte qu’elles reflètent du mieux possible le risque représenté par chaque assuré. On parlera ainsi de ‘segmentation tarifaire’. La segmentation se fait a priori, à l’entrée, mais aussi a posteriori, lors d’un renouvellement (on reviendra sur ce point). Mais au delà de l’élaboration d’une structure tarifaire, la segmentation peut être vu également comme le tracé d’une frontière entre ce que la compagnie souhaite et ne souhaite pas couvrir. La politique d’acceptation des risques ne peut se penser sans structure tarifaire. Et c’est précisément parce que la segmentation peut conduire à des exclusion qu’il est important de s’interroger sur les liens entre solidarité et assurance.

Mais la segmentation est surtout vu dans un sens plus large, où l’assureur peut, en plus de proposer des primes différentes, proposer des gammes de contrats variées (basées sur des franchises spécifiques, obligation de prendre des mesures de précautions particulières, exclusion de garanties). L’idée communément admise est que le choix d’une formule va révéler à l’assureur une information, comme le note Bourguignon (1994).

Continue reading Assurance, solidarité ou individualisation

Niveau et puissance de test (paramétrique)

Lors du dernier cours de statistique, nous avons vu comme se construisent des tests, et comment on compare leur puissance. Considérons pour cela des observations \{x_1,\cdots,x_n\} tirées suivant une loi uniforme \mathcal{U}([0,\theta]). Et on se pose la question de savoir si l’hypothèse H_0:\theta\leq\theta_0 est rejetée, ou acceptée, avec \theta_0=1 par exemple. L’hypothèse alternative sera ici le complémentaire, autrement dit H_1:\theta>\theta_0.

Comme on l’a vu, la construction d’un test se fait en utilisant une statistique de test (ou plus simplement un estimateur du paramètre inconnu), avec \psi(\mathbf{x})=\mathbf{1}(\widehat{\theta}(\mathbf{x}))\leq c. On prend la décision D_1 (consistant à rejeter H_0) si \psi(\mathbf{x})=1, et la décision D_0 (consistant à ‘accepter’ H_0) si \psi(\mathbf{x})=0.

Le risque d’erreur de première espèce est basé sur \mathbb{E}[\psi(\mathbf{X})]. Si l’hypothèse H_0 est une hypothèse simple H_0:\theta=\theta_0, le risque d’erreur de première espèce est https://latex.codecogs.com/gif.latex?\alpha=\mathbb{E}[\psi(\mathbf{X})\vert%20X_1,\cdots,X_n\sim%20F_{\theta_0}], soit \alpha=\mathbb{P}[\widehat{\theta}(\mathbf{X})>c \vert X_1,\cdots,X_n\sim F_{\theta_0}] pour un seuil c. On va déterminer le seuil c de manière à ce que le risque de première espèce soit 5%.

Pour cela, il nous fait un estimateur \widehat{\theta} dont on connaîtrait la loi. Commençons par l’estimateur du maximum de vraisemblance. Rappelons que l’estimateur du maximum de vraisemblance est \widehat{\theta}=\max\{x_i\}. Nous avons vu – dans la section sur les statistiques d’ordre – que \mathbb{P}[\max\{X_i\}]=\left(\frac{x}{\theta}\right)^n

pour x\leq\theta. Aussi, pour un seuil c, \mathbb{P}[\widehat{\theta}>c]=1-\left(\frac{c}{\theta}\right)^n

pour c\leq\theta. Le risque de première espèce \alpha=\underset{\theta\leq\theta_0}{\max}\left\lbrace\mathbb{P}\big[\widehat{\theta}\geq c\vert X\sim F_{\theta}\big]\right\rbrace

soit \alpha=\underset{\theta\leq\theta_0}{\max}\left\lbrace1-\left(\frac{c}{\theta}\right)^n\right\rbrace

qui atteindra son maximum en \theta_0. Aussi \alpha=1-\left(\frac{c}{\theta_0}\right)^n

soit c=\theta_0\cdot[1-\alpha]^{1/n}. Les régions d’acceptation et de rejet sont les suivantes,

Continue reading Niveau et puissance de test (paramétrique)

Regard personnel sur le Big Data

Sur variances.eu, la thématique Big Data du mois d’octobre se poursuit avec aujourd’hui le regard (personnel) de Xavier Dupré,

Je pianote ce soir sur le clavier de mon ordinateur. Il est 23h et je suis au hackathon organisé par Data Science Game – une association créée par d’anciens ENSAE. Les participants se débattent avec des données actuarielles pour prédire le taux d’achats de police d’assurance. Les 500 Mo de données auraient été considérées comme “big” il y a cinq ans. Elles tiennent maintenant aisément sur un ordinateur portable et il ne faut que quelques secondes pour calculer une moyenne. Tous ont déjà participé à des compétitions Kaggle, utilisé plusieurs librairies de machines learning. Certains obtiennent des résultats une heure après le début malgré des variables de types très variés. On se doute que ce n’est pas la première compétition à laquelle ils participent. Je circule d’une équipe à l’autre, les écrans montrent du python, un peu de R, une ligne de commande ouverte dans un coin, un notebook dans un autre, des machines distantes qui tournent quelque part. Ca discute beaucoup au sein des équipes. Elles sont pour la plupart composées de personnes qui se complètent, le data manager, le programmeur, le machine learner. Au dîner, j’échange avec quelques équipes. L’un d’eux me parle d’un de ses amis qui apprend des réseaux profonds pour entraîner une machine à jouer au jeu vidéo Doom. Il utilise des techniques de deep reinforcement learning. Nous abordons les domaines dans lesquels ces technologies s’appliqueront. La médecine sort en premier. Les transports en second. Les participants, vingt équipes de quatre, ont été sélectionnées sur une première compétition. Même si je joue le rôle de mentor, peu me posent des questions techniques. Il y a peu de choses qu’ils ne sachent déjà. Tous connaissent la programmation et les statistiques. Après 24h, certains assemblages deviennent plutôt complexes, la prédiction est issue d’un patchwork de modèles. On apprécie le coup de crayon des peintres, certains participants ont un bon coup de clavier. Les équipes prennent tour à tour la place de numéro 1.

[à suivre…]

La nécessaire transformation des métiers analytiques

Sur variances.eu, on continue la thèmatique Big Data du mois d’octobre, avec le point de vue de Baptiste Beaume, de Covéa, sur les aspects métiers.

L’arrivée puis la pérennisation du Big Data dans le monde professionnel impacte pleinement le marché de l’emploi : 130 000 créations de postes prévues en France d’ici 2020 (cf Plan Big Data) ou encore une pénurie annoncée, d’ici 2 ans, de 190 000 Data Scientist aux Etats Unis selon cf McKinsey. Ces profils à mi-chemin entre les statistiques et l’informatique devront être armés pour maîtriser les nouveaux types de données, les derniers outils ou langages mais aussi les méthodes pour exploiter toutes ces informations.

Ainsi de plus en plus de formations s’adaptent, voire se créent, pour proposer sur le marché ces CVs encore trop rares : ENSAE, Telecom Paris Tech, ENSAI, ENSIMAG, Polytechnique, Paris Dauphine, Telecom Nancy etc… Les étudiants d’aujourd’hui peuvent décider de prendre cette voie et ainsi s’assurer un bel avenir en termes d’employabilité – attention quand même à vivre une réelle histoire d’amour avec la donnée sans quoi le temps passé derrière un écran pourra sembler très long, la tension de ces profils sur le marché n’incitera pas les entreprises à les faire bifurquer vers d’autres domaines.

Les entreprises font actuellement toutes face à d’immenses enjeux de compétitivité liés à la maîtrise des données, la pénurie de Data Scientist n’en est que plus préoccupante dans une optique de recrutement. Cette problématique-là a toutefois tendance à en cacher une autre : quid des centaines de milliers de personnes qui ont déjà un métier analytique ? Des changements sont-ils à opérer ? Alors si oui lesquels ?

[à suivre….]

Unbiased Estimators vs. Minimizing a Quadratic Loss Function

Unbiased estimators are important in statistics. I guess because of Cramér Rao bound, for the variance. In the sense that if  , then  https://latex.codecogs.com/gif.latex?\text{Var}[\widehat{\theta}]\geq%20I_\theta^{-1}, where   denotes Fisher information (the proof was writen in an old post).

But what could we be the variance if   is not unbiased ?

Consider the following simple case, with a Gaussian i.id. sample   from a  . We know that the estimator of the Method of Moments is the same as the Maximum Likelihood estimator, i.e.  . And this estimator is efficient, in the sense that its variance is equal to Cramér-Rao lower bound,  .

But what if we consider another estimator? For instance  , with   not necessarily equal to 1. In that case, this estimator is (usually) biased since

 

And its variance is

 

We can visualise those two functions (the bias and the variance) using


n=10
alpha=seq(0,2,by=.01)
b=1-alpha
v=alpha^2/n
plot(alpha,b,xlab="alpha",col="red",type="l")
par(new=TRUE)
plot(alpha,v,col="blue",type="l",axes=FALSE
axis(4,)
mtext("bias",side=2,line=-1,col="red")
mtext("variance",side=4,line=-1,col="blue")

Observe that if is small (smaller than 1), the variance is smaller than the Cramér-Rao lower bound. And here, the mean squared error, defined as

which is, here,

The optimal value is obtained when the first order condition is satisfied

i.e.

So biased estimators can be more interesting than unbiased estimators, if the goal is the minimize the mean square error.

Re-parametrization and Maximum Likelihood

The maximum likelihood estimator is invariant in the sense that for all bijective function  , if   is the maximum likelihood estimator of   then  . Let  , then   is equal to  , and the likelihood function in   is  . And since   is the maximum likelihood estimator of  ,

hence,   is the maximum likelihood estimator of  .

For instance, the Bernoulli distribution is   with   and

 

Given sample  , the likelihood is

 

The log-likelihood is then

 

with ICI

 https://latex.codecogs.com/gif.latex?\frac{\partial}{\partial%20p}\log\mathcal{L}(p)=\frac{\sum%20x_i}{p}-\frac{n-\sum%20x_i}{1-p}.

Thus, the first order condition

 https://latex.codecogs.com/gif.latex?\frac{\partial}{\partial%20p}\log\mathcal{L}(p)=0

is satisfied when  . In order to illustrate, consider the following data


> set.seed(1)
> X=sample(0:1,size=15,replace=TRUE)
> X
[1] 0 0 1 1 0 1 1 1 1 0 0 0 1 0 1

The (negative) log-likelihood is here


> loglik=function(p){
+ -sum(log(dbinom(X,size=1,prob=p)))
+ }

that we can visualize below


> u=seq(0,1,by=.025)
> v=-Vectorize(loglik)(u)
> plot(u,v,type="l",xlab="",ylab="")

From calculations above, we know that the maximum likelihood estimator for  is


> mean(X)
[1] 0.5333333

The numerical version is


> (opt=optim(.5,loglik))
$par
[1] 0.5333008

$value
[1] 10.36385

$counts
function gradient
20 NA

$convergence
[1] 0

$message
NULL

Somehow, we were lucky here, because we did not say that the optimization was on the interval . Nevertheless, our estimator for the probability belongs to . In order to insure that the optimal value is in , we can consider some constrained optimization routine


> constrOptim(.5, loglik, grad=NULL,ui=matrix(c(1,-1),2,1), ci=c(0,-1))
$par
[1] 0.5333008

$value
[1] 10.36385

$counts
function gradient
20 NA

$convergence
[1] 0

$message
NULL

$outer.iterations
[1] 2

$barrier.value
[1] 6.909277e-05

On the previous graph, we did – indeed – reach that maximum of the log-likelihood


> abline(v=opt$par,col="red")

An alternative is to consider   (as in the exponential family). The log-likelihood is then

 

since

 

Here

 

Thus, the first order condition

 

is satisfied when

i.e.

 

From a numerical perspective, we have the same optimal value


> loglik=function(theta){
+ -sum(log(dbinom(X,size=1,prob=exp(theta)/(1+exp(theta)))))
+ }
> (opt=optim(0,loglik))
$par
[1] 0.1335938

$value
[1] 10.36385

$counts
function gradient
20 NA

$convergence
[1] 0

$message
NULL
> exp(opt$par)/(1+exp(opt$par))
[1] 0.5333489

Big Data (3) – Comment garantir au mieux la confidentialité ?

Dernier billet de la série publiée par Philippe Tassi, en ligne sur variance.eu, avec aujourd’hui une réflexion sur la confidentialité.

Les données ou statistiques, détenues ou élaborées par des administrations ou des entreprises, sont en général construites à partir d’informations individuelles, ce qui pose la question de la protection des sources, c’est-à-dire de la vie privée, compte tenu des progrès constants de la science et des process de traitement. Comment établir et maintenir la confiance du grand public, partie prenante numéro 1, tout en respectant l’équilibre entre promesse de confidentialité et utilisation des données recueillies ?

Pour y répondre, deux approches sont complémentaires : l’une est réglementaire, et montre que les Etats ont pris conscience depuis longtemps de la nécessité d’établir des garde-fous juridiques ; l’autre vise à s’appuyer sur la technologie en mettant en place des obstacles techniques empêchant la diffusion de données contre le gré de leur auteur.

[à suivre…]

Big Data (2) – phénomène éphémère ou de fond ?

Suite de la réflexion de Philippe Tassi sur variance.eu, avec aujourd’hui une réflexion sur phénomène ephémère ou de fond ?

Il est légitime de se poser cette question, tant certains sujets ont parfois été montés en épingle, avant de reperdre leur intérêt. Cependant, il ne semble pas que les données massives soient une mode passagère. Les vraies questions pourraient porter sur sa durée et sur ses conséquences sociétales. Combien d’années va vivre ce nouveau paradigme de juxtaposition de sources ? Comment se définir dans une société d’algorithmes traitant en temps réel des images numériques avec le minimum d’interventions humaines ?

Il est probable que nous soyons entrés pour au moins une décennie dans un monde à dominante scientifique. Non seulement les investissements sont réalisés dans la durée par les entreprises, et ne sont pas éphémères, mais s’y rajoute une réelle orientation politique.

En introduction, nous considérions que le numérique a donné encore plus de poids aux méthodologies, aux modélisations et aux technologies. Cela ne semble pas prêt de cesser. Non seulement parce que de nombreux secteurs d’activité développent ce type de sources – en 2014, par minute, le monde a créé 1,7 million de milliards d’octets de données, c’est-à-dire l’équivalent du contenu de 360 000 DVD –, mais aussi parce que c’est un espoir de croissance économique et d’emplois.

[à suivre…]

Les dérives du principe de précaution

« Dans le doute, abstiens-toi » dit la sagesse populaire. Le principe de précaution (en allemand Vorsorgeprinzip) est né de l’idée qu’il convient d’accepter qu’il existe un doute, ou une incertitude (scientifique), dans la connaissance des risques. Il y a un peu plus de 20 ans, la loi Barnier introduisait le principe de précaution dans le droit français, pour le « risque de dommages graves et irréversibles à l’environnement », communément appelés “risque environnemental”. Il y a un peu plus de 10 ans, il a été inscrit dans la Constitution, approuvé par 531 députés, exprimant ainsi un très large consensus politique, probablement aussi social. Mais aujourd’hui, le principe de précaution est évoqué dans des contextes tout aussi divers que le risque d’actes de terrorisme, mais aussi des procédures de droit civil, ou pénal. Quelles sont les conséquences de cette dérive de l’utilisation du principe de précaution ?

Continue reading Les dérives du principe de précaution

Big Data (1) – la fin des échantillons ou nouveau paradigme ?

Toujours dans la série de billets autour du big data sur le site variances.eu, le premier article d’une série de trois, par Philippe Tassi, sur le thème la fin des échantillons ou nouveau paradigme.

Depuis les premières tablettes cunéiformes du pays de Sumer jusqu’à la fin du XIXème siècle, on peut considérer que le monde a vécu sous le règne quasi-exclusif de « l’exhaustivité ».

Il y eut bien, au milieu du XVIIème siècle, de premières approches d’échantillonnage avec l’école dite de l’arithmétique politique de John Graunt et William Petty en Angleterre, et les avancées de Sébastien Le Prestre de Vauban et Jean-Baptiste Colbert en France : Vauban publie en 1686 un livre intitulé « Méthode généralle et facille pour faire le dénombrement des peuples ».

Cependant, ces tentatives innovantes ne se sont pas généralisées. La norme reste la vision exhaustive ; les pays les plus développés sur le plan statistique – n’oublions pas que statistique vient du latin status, mot dont l’un des sens est « attitude, posture, position, état » mais aussi qui possède un acception plus large en tant que « gouvernement, Etat » – qui gèrent des recensements, qu’ils soient de population, d’entreprises ou d’établissements industriels, ou de terres agricoles.

Le XXème siècle est marqué par le lent recul de l’exhaustif et la montée de plus en plus affirmée du paradigme de l’échantillonnage. On peut considérer que l’acte fondateur en est la communication d’Anders N. Kiaer, directeur du Bureau Central de Statistique du Royaume de Norvège, lors du Congrès de Berne de l’Institut International de Statistique d’août 1895 (IIS, maintenant dénommée ISI dans sa version anglaise International Statistical Institute). Intitulée Observations et expériences concernant des dénombrements représentatifs, elle ouvre un premier et long débat sur l’échantillonnage, la « pars pro toto ».

Au niveau le plus général, il faut apprécier tout particulièrement les quatre premières lignes qui composent le préambule de l’exposé de A. N. Kiaer, pertinentes et visionnaires.

[à suivre…]

Maximum de Vraisemblance, approche numérique

Histoire de voir des choses plus amusantes dans le cours de statistique, on va aborder le cas où le maximum de la vraisemblance n’a pas d’expression analytique simple, pour insister sur les aspects numériques. Par exemple, si on a un mélange de lois, comme un mélange de lois normales,

Cette loi peut être utile, comme loi paramétrique, pour décrire la loi de la taille de gens. En guise d’illustration, on pourra utiliser les données déjà utilisée au début du cours

> Davis=read.table(
+ "http://socserv.socsci.mcmaster.ca/jfox/Books/Applied-Regression-2E/datasets/Davis.txt")
> Davis[12,c(2,3)]=Davis[12,c(3,2)]
> X=Davis$height

La première étape est d’écrire la fonction de log-vraisemblance du mélange (ou plutôt l’opposé de la log-vraisemblance, les fonctions d’optimisation de R cherchant des minimums de fonctions)

> logL=function(theta){
+ p=theta[1]
+ m1=theta[2]
+ s1=theta[3]
+ m2=theta[4]
+ s2=theta[5]
+ logL=-sum(log(p*dnorm(X,m1,s1)+(1-p)*dnorm(X,m2,s2)))
+ return(logL)
+ }

La fonction la plus simple, pour maximiser la log-vraisemblance est la suivante (à partir d’un jeu de paramètre initiaux, histoire d’avoir un point de départ pour une descente de gradient)

> opt<-optim(c(.5,160,1,180,1),logL) > (theta=opt$par)
[1] 0.5987635 165.2547700 5.9410993 178.4856961 6.3547038

On reviendra dessus en cours car on peut faire ‘mieux’ en utilisant un algorithme d’optimisation sous contrainte, pour tenir compte du fait, par exemple, que la probabilité est forcément comprise entre 0 et 1 (c’est le cas ici, mais on a eu de la chance….).

Pour visualiser la densité estimée, on utilise

> hist(X,col="light green",probability=TRUE)
> lines(density(X))
> p=theta[1]
> m1=theta[2]
> s1=theta[3]
> m2=theta[4]
> s2=theta[5]
> v=p*dnorm(u,m1,s1)+(1-p)*dnorm(u,m2,s2)
> lines(u,v,col="blue",lwd=2)

Une autre solution est d’utiliser l’algorithme EM (sur lequel on reviendra en cours). Le point de départ part du constat que le problème de notre mélange est qu’on ne sait pas à quelle classe appartiennent nos gens. On va partir de valeurs initiales pour nos paramètres, et comparer les chances d’appartenir à chacunes des classes

> theta=c(.5,150,1,200,1)
> p1=dnorm(X,theta[2],theta[3])
> p2=dnorm(X,theta[4],theta[5])
> p=p1/(p1+p2)

A partir de ces probabilités d’appartenir à chacune des classes, on va estimer les paramètres de nos deux lois normales. Soit à l’aide de la formule fermée, soit en faisant du maximum de vraisemblance

> m1=sum(p*X)/sum(p)
> m2=sum((1-p)*X)/sum(1-p)
> s1=sqrt(sum(p*(X-m1)^2)/sum(p))
> s2=sqrt(sum((1-p)*(X-m2)^2)/sum(1-p))
> MlogL=function(a){
+ m1=a[1]
+ s1=a[2]
+ m2=a[3]
+ s2=a[4]
+ logL=-sum(log(p*dnorm(X,m1,s1)+(1-p)*dnorm(X,m2,s2)))
+ return(logL)
+ }
> opt<-optim(theta[2:5],MlogL)

L’idée est en fait d’avoir une boucle: on estime les probabilités d’appartenir dans les classes (compte tenu des paramètres des lois normales), et une fois qu’on a ces probabilités, on peut réestimer les paramètres des lois. Puis on recommence,

> theta=c(.5,150,1,200,1)
> for(s in 1:100){
+ p1=dnorm(X,theta[2],theta[3])
+ p2=dnorm(X,theta[4],theta[5])
+ p=p1/(p1+p2)
+ m1=sum(p*X)/sum(p)
+ m2=sum((1-p)*X)/sum(1-p)
+ s1=sqrt(sum(p*(X-m1)^2)/sum(p))
+ s2=sqrt(sum((1-p)*(X-m2)^2)/sum(1-p))
+ theta=c(mean(p),m1,s1,m2,s2)
+ }

On récupère alors les paramètres ‘optimaux’ de notre loi mélange

> hist(X,col="light green",probability=TRUE)
> lines(density(X))
> u=seq(130,210,by=.25)
> p=theta[1]
> m1=theta[2]
> s1=theta[3]
> m2=theta[4]
> s2=theta[5]
> v=p*dnorm(u,m1,s1)+(1-p)*dnorm(u,m2,s2)
> lines(u,v,col="red",lwd=2)

qui est relativement proche de notre estimation globale.

On reviendra au prochain cours sur ces méthodes ‘pratiques’ de calculs du maximum de vraisemblance.