Bootstrap tutorial in R [Eng]

Problème à signaler:

Télécharger



★★★★★★★★★★3.5 étoiles sur 5 basé sur 1 votes.
Votez ce document:

Bootstrap tutorial in R [Eng]

...

1 Intervalles de confiance bootstrap avec erreurs standard

Le manuel décrit comment construire un intervalle de confiance pour un paramètre de population en construisant un intervalle centré sur une estimation ponctuelle avec une marge d'erreur égale à deux fois l'erreur standard. Ici, nous allons estimer la taille de l’erreur type en appliquant le bootstrap et en échantillonnant de nombreux échantillons remplacés par l’échantillon d’origine, chacun ayant la même taille que l’échantillon d’origine, en calculant une estimation ponctuelle pour chacun et en recherchant l’écart type de cette distribution. des statistiques bootstrap.

"text-align: justify;">1.1 Atlanta Commute Times

L'ensemble de données CommuteAtlanta du manuel contient des variables relatives à un échantillon de 500 navetteurs dans la région d'Atlanta.

bibliothèque (Lock5Data)

data (CommuteAtlanta)

str (CommuteAtlanta)

## 'data.frame': 500 obs. de 5 variables:

## $ Ville: Facteur avec 1 niveau "Atlanta": 1 1 1 1 1 1 1 1 1 1 ...

## $ Age: int 19 55 48 45 48 43 48 41 47 39 ...

## $ Distance: int 10 45 12 4 15 33 15 4 25 1 ...

## $ Time: int 15 60 45 10 30 60 45 10 25 15 ...

## $ Sexe: Facteur avec 2 niveaux "F", "M": 2 2 2 1 1 2 2 1 2 1 ...

Ces navetteurs faisaient tous partie de l’enquête américaine sur le logement réalisée par le US Census Bureau et aucun ne travaillait de chez eux. Vous trouverez plus de détails sur la page d’aide du jeu de données.

? CommuteAtlanta

Pour construire l'intervalle de confiance pour le temps de déplacement moyen à Atlanta, nous devons trouver l'estimation ponctuelle (moyenne de l'échantillon) à partir de l'échantillon initial.

time.mean = avec (CommuteAtlanta, moyenne (Time))

time.mean

## [1] 29.11

Pour trouver l'erreur type, nous allons créer une énorme matrice avec 100

27 29 31

temps moyen

densité

Nous voyons une distribution pas trop asymétrique et en forme de cloche, plus ou moins. L'écart type de cette distribution est le suivant.

time.se = sd (boot.statistics)

time.se

## [1] 0,9414

Enfin, construisez l'intervalle de confiance. Ici, j'arrondis la marge d'erreur vers le haut et jusqu'à la première décimale pour qu'elle ait deux chiffres significatifs, et je fais preuve de prudence lorsque j'arrondis pour ne pas réduire l'intervalle trop petit.

me = plafond (10 * 2 * time.se) / 10

round (time.mean, 1) + c (-1, 1) * me

## [1] 27,2 31,0

Maintenant interprète en contexte.

Nous sommes à 95% sûrs que le temps de déplacement moyen à Atlanta parmi les navetteurs qui ne travaillent pas à domicile varie de 27,2 à 31 minutes.

L'exactitude de cette inférence dépend de la représentativité de l'échantillon initial par rapport à la population d'intérêt. Comme le US Census Bureau a probablement utilisé des méthodes d'échantillonnage précises, cela semble raisonnable. Nous pourrions examiner les références pour mieux comprendre les conditions associées à la procédure d’échantillonnage (en quelle année, en quelle saison, quel libellé a été utilisé pour collecter les données).

1.2 Ecrire une fonction

Comme il y a plusieurs étapes compliquées, il serait utile d'avoir une fonction pour toutes les faire afin qu'à l'avenir, nous puissions utiliser la fonction puis l'appeler. Voici un exemple de fonction qui prend un argument x supposé être un échantillon numérique et effectue le bootstrap B fois.

La fonction imprimera des informations utiles sur la console, dressera un graphique de la distribution de bootstrap et renverra les statistiques de bootstrap, l’intervalle, l’erreur standard et le graphique, sous forme de liste.

# Bret Larget

# 10 janvier 2014

# Une fonction de démarrage rapide pour un intervalle de confiance pour la moyenne

# x est un échantillon quantitatif unique

# B est le nombre souhaité d'échantillons de bootstrap à prendre

# binwidth est passé à geom_histogram ()

boot.mean = fonction (x, B, binwidth = NULL) {

n = longueur (x)

boot.samples = matrix (échantillon (x, taille = n * B, remplacer = VRAI), B, n)

boot.statistics = apply (boot.samples, 1, mean)

se = sd (statistiques de démarrage)

require (ggplot2)

if (is.null (binwidth))

binwidth = diff (range (boot.statistics)) / 30

p = ggplot (data.frame (x = boot.statistics), aes (x = x)) +

geom_histogram (aes (y = .. densité ..), binwidth = binwidth) + geom_density (color = "red")

parcelle (p)

intervalle = moyenne (x) + c (-1,1) * 2 * se

imprimer (intervalle)

return (list (boot.statistics = boot.statistics, intervalle = intervalle, se = se, tracé = p))

}

Voici comment utiliser la fonction pour les hauteurs des étudiants de notre cours.

étudiants = read.csv ("students.csv")

out = with (students, boot.mean (Height, B = 1000))

66 67 68 69 70

densité

## [1] 66,90 69,56

out $ interval

## [1] 66,90 69,56

Cependant, avant d'utiliser cet intervalle, nous devons faire preuve de prudence avant de conclure que la taille moyenne des étudiants de moins de 12 ans se situe entre 66,9 et 69,6 pouces, car l'échantillon d'élèves n'est pas aléatoire, mais constitue un exemple de commodité de notre classe. Voici deux variables de confusion possibles: le sexe et le pays d'origine. Le sexe d'un élève est-il associé à la taille et à la décision de prendre

Statistiques 302? Qu'en est-il du pays d'origine?

1.3 pour les boucles

Le guide R des auteurs implémente le bottillon en utilisant une boucle for. Plutôt que de prendre tous les échantillons à la fois, la boucle for en prend un à la fois. Généralement, le code R qui utilise apply () est plus efficace que le code qui utilise les boucles for. Essayez les deux pour un grand nombre de répliques bootstrap!

n = longueur (étudiants $ Hauteur)

B = 1000

résultat = rep (NA, B)

pour (i en 1: B) {

boot.sample = sample (n, replace = TRUE)

résultat [i] = moyenne (étudiants $ Height [boot.sample])

}

avec (étudiants, moyenne (taille) + c (-1, 1) * 2 * ds (résultat))

## [1] 66,89 69,58

1.4 Proportions

Considérez le problème de l’estimation de la proportion de pièces orange de Reese. J'ai choisi l'un des échantillons d'élèves au hasard (vraiment!) Et j'ai choisi un élève avec 11 bonbons orange et 19 bonbons sans orange.

Utilisons le bootstrap pour trouver un intervalle de confiance de 95% pour la proportion de pièces de Reese orange.

La solution la plus simple consiste à représenter les données de l’échantillon sous forme de vecteur avec 11 1 et 19 0 et à utiliser le même mécanisme que précédemment avec la moyenne de l’échantillon.

reeses = c (rep (1, 11), rep (0, 19))

reeses.boot = boot.mean (reeses, 1000, binwidth = 1/30)

0,2 0,4 0,6

densité

## [1] 0,1947 0,5386

Ainsi, sur la base de cet échantillon unique, nous sommes convaincus à 95% que la proportion réelle d’orange

Bootstrap tutorial in R [Eng]

Les pièces de Reese se situent entre 0,19 et 0,54. Si nous avions combiné les 48 échantillons en un seul, nous aurions pu résoudre le problème. Il y avait un total de 741 bonbons à l'orange et 699 ceux non oranges pour une proportion observée de 0,515.

reeses = c (rep (1, 741), rep (0, 699))

reeses.boot = boot.mean (reeses, 1000, binwidth = 0.005)

densité

## [1] 0,4888 0,5404

Avec cette taille d'échantillon beaucoup plus grande, le SE est beaucoup plus petit et l'intervalle de confiance est plus étroit. Quelqu'un d'autre pense que la vraie proportion est 0.5?

1.5 Différences de moyens

Je vais utiliser le jeu de données StudentSurvey du manuel pour illustrer l’utilisation du bootstrap pour estimer les différences de moyennes. Une variable intéressante est l'exercice, le nombre d'heures par semaine que chaque élève exerce.

data (StudentSurvey)

avec (StudentSurvey, résumé (exercice))

## Min. 1er Qu. Médiane Moyenne 3ème Qu. Max. NA

## 0.00 5.00 8.00 9.05 12.00 40.00 1

avec (StudentSurvey, résumé (genre))

## F M

## 169 193

avec (StudentSurvey, by (Exercise, Gender, mean, na.rm = TRUE))

## Sexe: F

## [1] 8.11

## ---------------------------------------------------- --------

## Sexe: M

## [1] 9.876

Nous voyons dans ce résumé que dans l’échantillon, les hommes exercent plus d’heures par semaine que les femmes.

Si nous traitons cet échantillon d'étudiants comme choisis au hasard parmi une population d'étudiants universitaires, nous pouvons estimer la différence de temps passé à faire de l'exercice pour chaque sexe. Notez que sans plus d'informations sur le processus d'échantillonnage, une telle inférence pourrait être faussée si les étudiants de l'échantillon diffèrent considérablement de la population.

Un élève n'a pas signalé de valeur pour l'exercice. Nous allons créer un nouveau bloc de données à utiliser qui élimine cet individu. Notez l'utilisation des crochets pour prendre un sous-ensemble des lignes (et de toutes les colonnes), ainsi que de! Is.na () pour rechercher tous les cas où la variable Exercise n'était pas manquante (NA).

newStudent = with (StudentSurvey, StudentSurvey [! is.na (Exercice),])

Avant de construire l’intervalle de confiance, voici un graphique des deux distributions. J'utilise des boîtes à moustaches côte à côte avec les points de données réels tracés dans un calque situé au-dessus, mais instables pour éviter une surexploitation. J'ai utilisé coord_flip () pour faciliter la comparaison des distributions et parce qu'un tracé court et large convient mieux à la page. J'ai souligné la boîte à moustaches en rouge et utilisé le rouge pour les couleurs aberrantes afin de pouvoir les distinguer des points de données tracés.

ggplot (newStudent, aes (x = sexe, y = exercice)) +

geom_boxplot (color = "red", outlier.colour = "red") +

geom_point (position = position_jitter (h = 0, w = 0.3)) +

ylab ('Nbre d'heures d'exercice par semaine') +

coord_flip ()

  • F

M

0 10 20 30 40

Nombre d'heures d'exercice par semaine

Le sexe

Nous utilisons length () pour trouver la taille d’échantillon de chaque groupe. Notez que le nombre de femelles est

n [1] et le nombre d'hommes est n [2].

n = with (newStudent, by (exercice, sexe, longueur))

n

## Sexe: F

## [1] 168

## ---------------------------------------------------- --------

## Sexe: M

## [1] 193

Le bloc de code suivant crée une matrice pour les hommes et une autre pour les femmes des échantillons avec remplacement de la même taille pour chacun. Nous utilisons ensuite apply () pour affiner la moyenne de chaque échantillon et utiliser les différences (hommes moins femmes) pour obtenir la distribution des statistiques bootstrap. Nous représentons ce graphique pour vérifier la symétrie et la forme d'une cloche. Tout a l'air bien.

B = 1000

female.samples = with (newStudent, matrix (sample (Exercice [Genre == "F")], taille = n [1] *

B, remplace = VRAI), B, n [1]))

male.samples = with (newStudent, matrix (sample (Exercice [Genre == "M")], taille = n [2] *

B, remplace = VRAI), B, n [2]))

female.means = appliquer (échantillons féminins, 1, moyenne)

male.means = appliquer (échantillons masculins, 1, moyenne)

boot.stat = male.means - female.means

ggplot (data.frame (x = boot.stat), aes (x = x)) + geom_density ()

Pour finir, prenez l’estimation ponctuelle (la différence dans les moyennes des échantillons) et ajoutez et soustrayez deux fois l’erreur type. Après avoir examiné la version non arrondie, arrondissez à une décimale pour deux chiffres significatifs.

xbars = with (newStudent, by (Exercise, Gender, mean))

moi = 2 * sd (boot.stat)

(xbars [2] - xbars [1]) + c (-1, 1) * me

## [1] 0,5397 2,9914

round ((xbars [2] - xbars [1]) + c (-1, 1) * moi, 1)

## [1] 0.5 3.0

1.6 Le paquet de démarrage

Il existe un package de démarrage avec une fonction boot () qui effectue le bootstrap dans de nombreuses situations. Je reviendrai sur l'exemple du Atlanta Commute Times. La fonction boot () nécessite trois arguments:

(1) les données de l'échantillon d'origine (une trame de données, une matrice ou un tableau); (2) une fonction pour calculer les statistiques à partir des données, le premier argument étant les données et le deuxième argument étant les indices des observations de l'échantillon bootstrap; (3) le nombre de réplicats bootstrap.

Un exemple vous aidera à comprendre. bibliothèque (boot)

data (CommuteAtlanta)

my.mean = function (x, indices) {

retour (moyenne (x [indices]))

}

time.boot = boot (CommuteAtlanta $ Time, my.mean, 10000)

Notez que my.mean (CommuteAtlanta $ Time, 1: length (CommuteAtlanta $ Time) calcule la moyenne de l'échantillon d'origine.

L'objet time.boot est une liste avec de nombreux éléments. L'un est time.boot $ t0 qui est la moyenne d'échantillon des données d'origine. Time.boot $ t est une autre collection de statistiques d’amorçage utilisables comme ci-dessus. Mais la fonction intégrée boot.ci () calculera les intervalles de confiance au démarrage à l’aide de plusieurs méthodes.

boot.ci (time.boot)

## CALCULS D'INTERVALLE DE CONFIANCE DE BOOTSTRAP

## Basé sur 10000 réplicats d'amorçage

## APPEL :

## boot.ci (boot.out = time.boot)

## intervalles:

## Niveau Normal Basique

## 95% (27,28, 30,96) (27,24, 30,92)

## Niveau Percentile BCa

## 95% (27.30, 30.98) (27.40, 31.13)

## Calculs et intervalles sur l'échelle originale

Basic utilise l'erreur standard estimée. Le centile utilise les centiles. La BCa utilise également des centiles, mais elle a été ajustée pour tenir compte des biais et des biais.

2 Percentile Bootstrap à venir. . .

L'idée des intervalles de confiance utilisant les centiles du bootstrap consiste à sélectionner les points d'extrémité situés au milieu de la distribution du bootstrap correspondant au niveau de confiance souhaité. Voici un exemple utilisant l’exercice 3.109.

require (Lock5Data)

données (ImmuneTea)

thé = avec (ImmuneTea, InterferonGamma [Drink == "Thé"])

café = avec (ImmuneTea, InterferonGamma [Drink == "Café"]])

thé.moyen = moyenne (thé)

café.moyen = moyenne (café)

Bootstrap tutorial in R [Eng]

tea.n = longueur (thé)

café.n = longueur (café)

B = 100000

# crée des tableaux vides pour les moyennes de chaque échantillon

tea.boot = numérique (B)

coffee.boot = numérique (B)

# Utilisez une boucle for pour prélever les échantillons pour (i in 1: B) {

tea.boot [i] = moyenne (échantillon (thé, taille = thé.n, remplacer = VRAI))

coffee.boot [i] = moyenne (échantillon (café, taille = café.n, remplace = VRAI))

}

boot.stat = tea.boot - coffee.boot

# Recherchez les points de terminaison pour les intervalles de confiance bootstrap de 90%, 95% et 99% à l'aide de centiles.

quantile (boot.stat, c (0.05,0.95))

## 5% 95%

# 4.018 29.764

quantile (boot.stat, c (0.025,0.975))

## 2,5% 97,5%

## 1.491 32.045

quantile (boot.stat, c (0.005,0.995))

## 0,5% 99,5%

## -3,745 36,264


1