Qualité de l'ajustement pour les histogrammes 2D

19

J'ai deux ensembles de données représentant les paramètres des étoiles: un observé et un modélisé. Avec ces ensembles, je crée ce qu'on appelle un diagramme à deux couleurs (TCD). Un échantillon peut être vu ici:

histogrammes

A étant les données observées et B les données extraites du modèle (sans parler des lignes noires, les points représentent les données) Je n'ai qu'un seul diagramme A , mais je peux produire autant de diagrammes B différents que je veux, et ce dont j'ai besoin est de garder celui qui correspond le mieux A .

J'ai donc besoin d'un moyen fiable de vérifier la qualité de l'ajustement du diagramme B (modèle) au diagramme A (observé).

En ce moment, je crée un histogramme ou une grille 2D (c'est ce que j'appelle, peut-être qu'il a un nom plus approprié) pour chaque diagramme en regroupant les deux axes (100 cases pour chacun) Ensuite, je passe par chaque cellule de la grille et je trouve la différence absolue en nombre entre A et B pour cette cellule particulière. Après avoir passé par toutes les cellules, je résume les valeurs de chaque cellule et je termine donc avec un seul paramètre positif représentant la qualité de l' ajustement ( ) entre A et B . Le plus proche de zéro, meilleur est l'ajustement. Fondamentalement, voici à quoi ressemble ce paramètre:gF

gF=jej|unejej-bjej|; où est le nombre d'étoiles dans le diagramme A de cette cellule particulière (déterminée par ) et est le nombre de chambres .unejejjejbjej

Voici à quoi ressemblent ces différences de nombre dans chaque cellule dans la grille que je crée (notez que je n'utilise pas de valeurs absolues de dans cette image mais je les utilise pour calculer le paramètre ):(unejej-bjej)(unejej-bjej)gF

hess

Le problème est que j'ai été informé que cela pourrait ne pas être un bon estimateur, principalement parce que, à part dire que cet ajustement est meilleur que cet autre parce que le paramètre est plus bas , je ne peux vraiment rien dire de plus.


Important :

(merci @PeterEllis d'avoir soulevé cette question)

1- Les points en B ne sont pas liés l' un à l' un des points en A . C'est une chose importante à garder à l'esprit lors de la recherche du meilleur ajustement: le nombre de points dans A et B n'est pas nécessairement le même et le test de qualité de l'ajustement devrait également tenir compte de cet écart et essayer de le minimiser.

2- Le nombre de points dans chaque ensemble de données B (sortie du modèle) que j'essaye d'adapter à A n'est pas fixe.


J'ai vu le test du Chi-Squared utilisé dans certains cas:

je(Oje-Eje)2/Eje ; où est la fréquence observée (modèle) et est la fréquence attendue (observation).OjeEje

mais le problème est: que dois-je faire si est nul? Comme vous pouvez le voir dans l'image ci-dessus, si je crée une grille de ces diagrammes dans cette plage, il y aura beaucoup de cellules où est nul.EjeEje

De plus, j'ai lu que certaines personnes recommandent d'appliquer un test de Poisson de vraisemblance logarithmique dans des cas comme celui-ci où des histogrammes sont impliqués. Si cela est correct, j'apprécierais vraiment que quelqu'un puisse m'instruire sur la façon d'utiliser ce test dans ce cas particulier (rappelez-vous, ma connaissance des statistiques est assez faible, alors s'il vous plaît restez aussi simple que possible :)

Gabriel
la source
Les points de B ont-ils une relation un à un avec les points de A (par exemple, chacun est une étoile particulière) ou est-ce plus abstrait que cela?
Peter Ellis
Hi @PeterEllis, aucun point B ne sont pas liés l' un à l' un des points en A . En fait, c'est une autre chose importante à garder à l'esprit lors de la recherche du meilleur ajustement: le nombre de points dans A et B n'est pas nécessairement égal.
Gabriel
Salut - question intéressante, je vais essayer d'écrire une bonne réponse. Chaque version de B a-t-elle le même nombre de points, ou varient-elles aussi?
Peter Ellis
Ils varient également, seul le nombre de points en A reste constant. Vous n'avez AUCUNE idée de combien vous m'aideriez si vous m'aidiez à comprendre cela @PeterEllis.
Gabriel
Cette question ressemble beaucoup à ce sujet: stats.stackexchange.com/questions/71036/… Où j'ai fourni une réponse.
L Fischman

Réponses:

14

OK, j'ai largement révisé cette réponse. Je pense que plutôt que de regrouper vos données et de comparer les nombres dans chaque groupe, la suggestion que j'avais enterrée dans ma réponse originale d'ajuster une estimation de densité de noyau 2D et de les comparer est une bien meilleure idée. Encore mieux, il y a une fonction kde.test () dans le paquet ks de Tarn Duong pour R qui rend cela aussi simple que possible.

Consultez la documentation de kde.test pour plus de détails et les arguments que vous pouvez modifier. Mais en gros, il fait à peu près exactement ce que vous voulez. La valeur p qu'il renvoie est la probabilité de générer les deux ensembles de données que vous comparez sous l'hypothèse nulle qu'ils étaient générés à partir de la même distribution. Donc, plus la valeur de p est élevée, meilleur est l'adéquation entre A et B.Voir mon exemple ci-dessous où cela ramasse facilement que B1 et A sont différents, mais que B2 et A sont plausiblement les mêmes (c'est ainsi qu'ils ont été générés) .

# generate some data that at least looks a bit similar
generate <- function(n, displ=1, perturb=1){
    BV <- rnorm(n, 1*displ, 0.4*perturb)
    UB <- -2*displ + BV + exp(rnorm(n,0,.3*perturb))
    data.frame(BV, UB)
}
set.seed(100)
A <- generate(300)
B1 <- generate(500, 0.9, 1.2)
B2 <- generate(100, 1, 1)
AandB <- rbind(A,B1, B2)
AandB$type <- rep(c("A", "B1", "B2"), c(300,500,100))

# plot
p <- ggplot(AandB, aes(x=BV, y=UB)) + facet_grid(~type) + 
    geom_smooth() +     scale_y_reverse() + theme_grey(9)
win.graph(7,3)
p +geom_point(size=.7)

entrez la description de l'image ici

> library(ks)
> kde.test(x1=as.matrix(A), x2=as.matrix(B1))$pvalue
[1] 2.213532e-05
> kde.test(x1=as.matrix(A), x2=as.matrix(B2))$pvalue
[1] 0.5769637

MA RÉPONSE ORIGINALE CI-DESSOUS, EST CONSERVÉE UNIQUEMENT PARCE QU'IL Y A MAINTENANT DES LIENS D'AUTRE PART QUI NE FONT PAS SENS

Premièrement, il peut y avoir d'autres façons de procéder.

Justel et al ont proposé une extension multivariée du test de qualité d'ajustement de Kolmogorov-Smirnov qui, je pense, pourrait être utilisée dans votre cas, pour tester l'adéquation de chaque ensemble de données modélisées avec l'original. Je n'ai pas pu trouver une implémentation de ceci (par exemple dans R) mais peut-être que je n'ai pas regardé assez fort.

Alternativement, il peut y avoir un moyen de le faire en ajustant une copule à la fois aux données d'origine et à chaque ensemble de données modélisées, puis en comparant ces modèles. Il existe des implémentations de cette approche dans R et dans d'autres endroits, mais je ne les connais pas particulièrement, donc je n'ai pas essayé.

Mais pour répondre directement à votre question, l'approche que vous avez adoptée est raisonnable. Plusieurs points se suggèrent:

  • À moins que votre ensemble de données ne soit plus grand qu'il n'y paraît, je pense qu'une grille de 100 x 100 est trop de bacs. Intuitivement, je peux imaginer que vous concluez que les différents ensembles de données sont plus différents qu'ils ne le sont simplement en raison de la précision de vos bacs signifie que vous avez beaucoup de bacs avec un petit nombre de points, même lorsque la densité de données est élevée. Mais en fin de compte, c'est une question de jugement. Je vérifierais certainement vos résultats avec différentes approches du binning.

  • Une fois que vous avez effectué votre regroupement et que vous avez converti vos données en (en fait) une table de contingence de comptages avec deux colonnes et un nombre de lignes égal au nombre de bacs (10 000 dans votre cas), vous avez un problème standard de comparaison de deux colonnes des chefs d'accusation. Un test du chi carré ou l'ajustement d'une sorte de modèle de Poisson fonctionnerait, mais comme vous le dites, il y a une gêne en raison du grand nombre de dénombrements nuls. L'un ou l'autre de ces modèles est normalement ajusté en minimisant la somme des carrés de la différence, pondérée par l'inverse du nombre attendu de comptes; lorsque cela approche de zéro, cela peut causer des problèmes.

Edit - le reste de cette réponse, je ne pense plus que ce soit une approche appropriée.

ng×2

ng×2ng

J'ai simulé certaines données pour qu'elles ressemblent un peu aux vôtres et j'ai constaté que cette approche était assez efficace pour identifier les ensembles de données "B" qui avaient été générés à partir du même processus que "A" et lesquels étaient légèrement différents. Certainement plus efficace que l'oeil nu.

  • ng×2un problème si vous utilisez uniquement la somme des différences absolues ou des différences au carré, comme vous le proposiez à l'origine). Cependant, il importe que chacune de vos versions de B ait un nombre différent de points. Fondamentalement, de plus grands ensembles de données B auront tendance à renvoyer des valeurs de p plus faibles. Je peux penser à plusieurs solutions possibles à ce problème. 1. Vous pouvez réduire tous vos ensembles de données B à la même taille (la taille du plus petit de vos ensembles B), en prélevant un échantillon aléatoire de cette taille dans tous les ensembles B qui sont plus grands que cette taille. 2. Vous pouvez d'abord ajuster une estimation de densité de noyau bidimensionnelle à chacun de vos ensembles B, puis simuler des données de cette estimation de tailles égales. 3. vous pouvez utiliser une sorte de simulation pour déterminer la relation entre les valeurs de p et la taille et l'utiliser pour «corriger» les valeurs de p que vous obtenez à partir de la procédure ci-dessus afin qu'elles soient comparables. Il existe probablement aussi d'autres alternatives. Lequel dépendra de la façon dont les données B ont été générées, de la différence des tailles, etc.

J'espère que cela pourra aider.

Peter Ellis
la source
J'ai fait quelques corrections de faute mineures; J'espère que ça ne vous dérange pas. Il peut y avoir un moyen de formater les choses pour faire ressortir les idées principales un peu plus en évidence, en particulier au dernier point. Mais je ne voulais pas non plus être trop zélé. À votre santé. :)
Cardinal
pas de problème. J'ai eu du mal avec un bon moyen de formater mon dernier point - ce que je voulais, c'était une hiérarchie de liste numérotée sous un point. Mais je n'ai pas trouvé comment faire ça.
Peter Ellis
Ouais, j'ai aussi joué avec ça brièvement, parce que c'est ce que tu avais l'air de vouloir. Je ne pouvais pas comprendre rapidement comment faire cela et j'étais très hésitant à apporter des modifications en gros à la mise en page, alors j'ai pensé que je ferais simplement un commentaire à la place. :)
Cardinal
1
Je recommande le livre de Wilcox en tant que texte de statistiques générales qui utilise R (voir ma réponse stats.stackexchange.com/questions/25632/… ). Bien qu'il ne couvre pas le texte exact de Fisher, il y a suffisamment de détails sur le texte spécifique sur le Web qui auront plus de sens si vous avez des informations dans ce livre sur des tests similaires. Il pourrait y avoir un meilleur texte sur ce type spécifique de problème d'ajustement, mais je pense que le livre de Wilcox est génial en tant qu'introduction générale qui vous amène rapidement à un niveau élevé.
Peter Ellis
1
Sensationnel. Vous avez répondu au diable de cette chose. S'il y avait un «meilleur échange de pile», ce serait là.
Colin K