Obtention de la densité de données probables

Question

J'ai besoin d'analyser certaines données sur les sessions Internet pour une ligne DSL. Je voulais voir comment les durées de sessions sont distribuées. J'ai pensé qu'un moyen simple de le faire serait de commencer par faire un graphique de densité de probabilité de la durée de toutes les sessions. J'ai chargé les données dans R et j'ai utilisé la fonction density(). Donc, c'était quelque chose comme ça

plot(density(data$duration), type = "l", col = "blue", main = "Density Plot of Duration", 
    xlab = "duration(h)", ylab = "probability density") 

Je suis nouveau à R et ce genre d'analyse. C'est ce que j'ai trouvé en passant par google. J'ai eu une intrigue mais il me restait quelques questions. Est-ce la bonne fonction de faire ce que j'essaie de faire ou y a-t-il autre chose?

Dans le graphique, j'ai trouvé que l'échelle de l'axe Y était de 0 à 1,5. Je ne comprends pas comment ça peut être 1.5, ne devrait-il pas être de 0 ... 1?

En outre, je voudrais obtenir une courbe plus lisse. Depuis, l'ensemble de données est vraiment grand les lignes sont vraiment déchiquetées. Ce serait plus agréable de les aplanir quand je présente cela. Comment ferais-je cela?

Answer 1

Vous devez jouer avec le paramètre de bande passante (bw) pour modifier le lissage de la courbe. Généralement, R fait un bon travail et donne automatiquement une courbe agréable et lisse, mais ce n'est peut-être pas le cas pour votre ensemble de données spécifique.

En ce qui concerne l'appel que vous utilisez, oui, c'est correct, type="l" n'est pas nécessaire, c'est la valeur par défaut utilisée pour tracer les objets de densité. L'aire sous la courbe (c'est-à-dire l'intégrale de -Inf à + Inf de votre fonction de densité) sera = 1.

Maintenant, est-ce qu'une courbe de densité est la meilleure chose à utiliser dans votre cas? Peut-être, peut-être pas ... cela dépend vraiment du type d'analyse que vous voulez faire. Probablement en utilisant hist sera suffisant, et peut-être plus informatif que vous pouvez sélectionner des intervalles de durée spécifiques (voir ?hist pour plus d'informations).

Answer 2

Comme nico a dit, vous devriez vérifier hist, mais vous pouvez également combiner les deux. Ensuite, vous pouvez appeler la densité avec lines à la place. Exemple:

duration <- rpois(500, 10) # For duration data I assume Poisson distributed 
hist(duration, 
    probability = TRUE, # In stead of frequency 
    breaks = "FD",  # For more breaks than the default 
    col = "darkslategray4", border = "seashell3") 
lines(density(duration - 0.5), # Add the kernel density estimate (-.5 fix for the bins) 
    col = "firebrick2", lwd = 3) 

devrait vous donner quelque chose comme:

Notez que l'estimation de la densité du noyau suppose un noyau gaussien par défaut. Mais la bande passante est souvent le facteur le plus important. Si vous appelez directement density il signale la bande passante estimée par défaut:

> density(duration) 

Call: 
     density.default(x = duration) 

Data: duration (500 obs.);  Bandwidth 'bw' = 0.7752 

     x     y    
Min. : 0.6745 Min. :1.160e-05 
1st Qu.: 7.0872 1st Qu.:1.038e-03 
Median :13.5000 Median :1.932e-02 
Mean :13.5000 Mean :3.895e-02 
3rd Qu.:19.9128 3rd Qu.:7.521e-02 
Max. :26.3255 Max. :1.164e-01 

Ici, il est 0,7752. Vérifiez-le pour vos données et jouez avec lui comme suggéré par nico. Vous voudrez peut-être regarder ?bw.nrd.

Answer 3

J'allais ajouter ceci comme commentaire à la réponse précédente, mais c'est trop gros. Le biais apparent est dû à la manière dont les valeurs sont regroupées dans un histogramme. C'est souvent une erreur d'utiliser des histogrammes pour des données discrètes. Voir ci-dessous ...

set.seed(1001) 
tmpf <- function() { 
    duration <- rpois(500, 10) # For duration data I assume Poisson distributed 
    hist(duration, 
     probability = TRUE, # In stead of frequency 
     breaks = "FD",  # For more breaks than the default 
     col = "darkslategray4", border = "seashell3", 
     main="",ann=FALSE,axes=FALSE,xlim=c(0,25),ylim=c(0,0.15)) 
    box() 
    lines(density(duration), # Add the kernel density estimate 
     col = "firebrick2", lwd = 3) 
    par(new=TRUE) 
    plot(table(factor(duration,levels=0:25))/length(duration), 
     xlim=c(0,25),ylim=c(0,0.15),col=4,ann=FALSE,axes=FALSE) 
} 

par(mfrow=c(3,3),mar=rep(0,4)) 
replicate(9,tmpf())