- Importer la table depuis un entrepôt distant

MaTableEnLigne <- read.table("https://git.io/vbiaS" , header = TRUE)
# J'affecte les données entreposées sur un serveur distant à un objet R local que je nomme 'MaTableEnLigne'. 
# Les noms des objets R ne doivent pas commencer par un chiffre ni contenir des espaces ou des caractères spéciaux.
# 'header=TRUE' signifie que la 1ère ligne contient les noms des colonnes. Pour avoir la description
# d'une fonction R et de ses arguments, il faut faire:  help(nom de la fonction). Par exemple:  help(read.table)

- Écrire cette table au format texte (.txt) dans le dossier de travail

write.table(MaTableEnLigne , "MaTableEnLocal.txt" , row.names = FALSE)

- Importer la table au format texte depuis le dossier de travail

MaTableEnLocal <- read.table("MaTableEnLocal.txt" , header = TRUE) 

- Vérification utile

identical(summary(MaTableEnLocal) , summary(MaTableEnLigne))

## [1] TRUE

- Exporter une table depuis R au format Excel vers le dossier de travail local

write.xlsx(MaTableEnLigne , "MaTableEnLocal.xlsx")

- Importer dans R une table au format Excel depuis le dossier de travail local

CM1 <- read.xlsx("MaTableEnLocal.xlsx" , sheet = 1)

- Réaffecter à un objet existant une autre table

CM1 <- MaTableEnLigne # Je réaffecte à l'objet CM1 le contenu de la table distante.

colnames(CM1) # Noms des variables de la table.

## [1] "x" "y"

CM1$y # l'opérateur '$' entre 'CM1' et 'y' permet d'extraire de la table 'CM1' la colonne 'y'.

##  [1]  8.935373  6.523501 18.649757 22.613518 15.965711 16.377628 22.271611
##  [8] 19.426137 16.603166 25.587190 23.682441 26.848028 23.655190 30.597357
## [15] 32.777300 32.020974 27.205078 33.229944 35.586207 42.126671

sum(CM1$y > 20) # Nombre de valeurs de y supérieures à 20.

## [1] 13

sum(CM1$y > 20)/nrow(CM1) # La proportion des valeurs>20 dans y.

## [1] 0.65

plot(CM1$x , CM1$y , col = ifelse(CM1$x < 10 , 'magenta' , 'gold') , 
 pch = 19 , cex = 2 , main = 'Graphique avec R: couleur magenta pour les valeurs de x inf. à 10.')
# Repérer les individus sur le graphique.
text(CM1$x , CM1$y , rownames(CM1) , col ='white' , cex =0.7) #la fonction 'text()' superpose du texte sur le graph. 

suppressPackageStartupMessages({library(plotly)})
library(plotly)
plot_ly(CM1 , x = ~x , y = ~y , type='scatter' , mode='markers')%>%
layout(title = "Ce paragraphe n'est pas au programme , c'est une démonstration")

y observés: \(y_i = \hat{\beta}_0 + \hat{\beta_1} x_{i} + \hat{\epsilon}_i\)
y estimés : \(\hat{y}_i = \hat{\beta_0} + \hat{\beta_1} x_{i}\)

\(\hat{\beta}_0\) L’ordonnée à l’origine et le test de l’ordonnée à l’origine nulle.

\(\hat{\beta_1}\) La pente et le test de pente nulle.

\(\hat{\epsilon}\) Les résidus, leur écart-type et degrés de liberté.

\(R^2\) Le coefficient de détermination et \(R^2_{adjust}\) le coefficient de détermination ajusté.

\(\hat{\beta_0}\) et \(\hat{\beta_1}\) sont des estimateurs sans biais de \(\beta_0\) et \(\beta_1\) respectivement

Calculs intermédiaires : moyenne et variance de x, moyenne et variance de y et la cov(x,y).

\[\mu_x = \frac1n\sum_{i=1}^{n} x_i\]

mx <- mean(CM1$x) ; mx #Moyenne de x.

## [1] 9.239499

\[\mu_y = \frac1n\sum_{i=1}^{n} y_i\]

my <- mean(CM1$y) ; my #Moyenne de y.

## [1] 24.03414

n  <- nrow(CM1) ; n  # Nombre d'observations (Ici on n'a pas de valeurs manquantes).

## [1] 20

\[Cov(x,y) = \frac1n\sum_{i=1}^n(x_i-\mu_x)*(y_i-\mu_y)\]

cov_xy  <- sum((CM1$x-mx)*(CM1$y-my))/(n-1) ; cov_xy #Covariance de x et de y.

## [1] 39.63647

\[Var(x) = \frac1n\sum_{i=1}^n(x_i-\mu_x)^2 \]

var_x   <- sum((CM1$x-mx)^2)/(n-1) ; var_x #Variance de x.

## [1] 22.36218

\[ Var(y) = \frac1n\sum_{i=1}^n(y_i-\mu_y)^2 \]

var_y   <- sum((CM1$y-my)^2)/(n-1) ; var_y #Variance de y.

## [1] 79.90522

\[\hat{\beta_1} = \frac{Cov(x,y)}{Var(x)}\]

Beta1   <- cov_xy/var_x ; Beta1 #Beta 1.

## [1] 1.772478

\[\hat{\beta_0} = \mu_y - \hat{\beta_1} * \mu_x\]

Beta0   <- my-Beta1*mx ; Beta0 #Beta 0.

## [1] 7.65733

- Vérification graphique des Coefficients \(\hat{\beta_0}\) et \(\hat{\beta_1}\)

plot(CM1$x , CM1$y)
abline(a = Beta0 , b = Beta1)

- Les résidus

\[\hat{\epsilon}_i = y_i - \hat{\beta_0} + \hat{\beta_1} x_{i}\]

y_obs   <- CM1$y #y observé (les données réelles de y).

\[\hat{y}_i = \hat{\beta_0} + \hat{\beta_1} x_{i} \]

y_est   <- Beta0 + Beta1*CM1$x #y estimé par le modèle.

Vérification

sum(y_obs) == sum(y_est) # Si notre calcul de 'y_est' est correct alors le résultat doit être 'TRUE'.

\[\hat{\epsilon}_i = y_i - \hat{y}_i \quad \quad ; \quad \quad \sum_{i=1}^n \hat{\epsilon}_i = 0 \]

residus <- y_obs - y_est ; sum(residus) #calcul des résidus, leur somme est = 0.

## [1] 1.687539e-14

Décomposition de la variance : \(SC_{total} = SC_{expli} + SC_{resid}\)

\[SC_{total} = \sum_{i=1}^n(y_i-\mu_y)^2 \]

SC_total <- sum((y_obs - my)    ^2) ; SC_total #Somme des carrés totale.

## [1] 1518.199

\[ SC_{expli} = \sum_{i=1}^n(\hat{y}_i-\mu_y)^2 \]

SC_expli <- sum((y_est - my)    ^2) ; SC_expli #Somme des carrés due à la régression.

## [1] 1334.841

\[ SC_{resid} = \sum_{i=1}^n(y_i - \hat{y}_i)^2 \]

SC_resid <- sum((y_obs - y_est) ^2) ; SC_resid #Somme des carrés résiduelle.

## [1] 183.3586

- Degrés de liberté.

n = : nombre d’individus = taille de l’échantillon (sans les valeurs manquantes) et r = Nombre de coefficients estimés (intercept et pente). Donc r = 2.

n = nrow(CM1) # Vrai s'il n'y a pas de valeurs manquantes dans les 2 variables du modèle.
r = 2

ddl_SC_total <-  n - 1              ; ddl_SC_total    # n - 1

## [1] 19

ddl_SC_resid <-  n - r              ; ddl_SC_resid    # n - r 

## [1] 18

ddl_SC_expli <- (n - 1) - (n - r)   ; ddl_SC_expli    # r - 1 

## [1] 1

- écart-type des résidus (Residual standard error)

\(\hat{\sigma}^2 = \frac{\sum_{i=1}^n\hat{\epsilon}_i^2}{n-2}\) est un estimateur sans biais de \(\sigma^2\)

sigma2 <- SC_resid/ddl_SC_resid #variance (résidus)
SC_x   <- sum((CM1$x - mean(CM1$x))^2) # somme des carrées de x

Residual standard error

sqrt(sigma2)

## [1] 3.191644

Estimate Std. Error Beta0 \[ Var(\hat{\beta}_0) = \sigma^2(\frac1n+\frac{\mu_x^2}{\sum_{i=1}^n(x_i-\mu_x)^2}) \]

B0_std_err <-  sqrt(sigma2*(1/n+mean(CM1$x)^2/sum((CM1$x-mean(CM1$x))^2))) ; B0_std_err

## [1] 1.598764

P Valeur : Sous \(H_0\) \[ Z_0 = \frac{\hat{\beta}_0 -\beta_0}{Std.Error(\hat{\beta}_0)} \sim Student(n-2) \]

Z0 <- Beta0/B0_std_err # Z0 : la statistique de test de l'ordonnée à l'origine nulle. Test bilatéral.
p_value_Z0 <- 2*(pt(Z0 , n - 2 , lower.tail = F)); p_value_Z0

## [1] 0.0001467133

Estimate Std. Error Beta1

\[Var(\hat{\beta}_1) = \frac{\sigma^2}{\sum_{i=1}^{n}(x_i-\mu_x)^2}\]

B1_std_err <-  sqrt(sigma2/SC_x) ; B1_std_err

## [1] 0.1548391

P Valeur : Sous \(H_0\) \[ Z_1 = \frac{\hat{\beta}_1 -\beta_1}{Std.Error(\hat{\beta}_1)} \sim Student(n-2) \]

Z1 <- Beta1/B1_std_err # Z1 : Calcul de la p valeur du test de pente nulle. Test bilatéral.
p_value_Z1 <- 2*(pt(Z1 , n - 2 , lower.tail = F)); p_value_Z1

## [1] 1.074088e-09

\[ R^2 = \frac{SC_{expli}}{SC_{total}} =\frac{\sum_{i=1}^n(\hat{y}_i-\mu_y)^2}{\sum_{i=1}^n(y_i-\mu_y)^2} \]

R2 <- SC_expli/SC_total ; R2 # coefficient de détermination (dans l'échantillon)

## [1] 0.8792262

R2_ajust <- 1 -((n-1)/(n-r))*(1-R2) ; R2_ajust # # R2 ajusté (dans la population)

## [1] 0.8725166

\[ \operatorname{Cor}(x,y) = \frac{\operatorname{Cov}(x,y)}{\sigma_x \sigma_y} \]

R <- cov_xy/sqrt(var_x*var_y) ; R == sqrt(R2) # Coefficien de corrélation

## [1] TRUE

- Test de Fisher de comparaison de deux variances F-statistic:

\[F=\frac{SC_{expli}/1}{SC_{resid}/(n-2)} \sim Fisher(1,n-2) \]

\[F=\frac{SC_{expli}/ddl_{expli}}{SC_{resid}/ddl_{resid}} \]

Fisher1  <- (SC_expli/ddl_SC_expli)/(SC_resid/ddl_SC_resid) ; Fisher1

## [1] 131.039

Lien entre le R2 et la statistique de Fisher

\[F = \frac{ddl_{SC.res}}{1} \frac{R^2}{1-R^2}\]

Fisher2  <- (ddl_SC_resid*R2) / (ddl_SC_expli*(1-R2)) ;Fisher2

## [1] 131.039

## P valeur
pf(Fisher1, ddl_SC_expli , ddl_SC_resid , lower.tail = F)

## [1] 1.074088e-09

antoine <- lm(CM1$y ~ CM1$x) #J'appelle le modèle créé 'antoine', y est la réponse, x est le prédicteur.
summary(antoine)  

## 
## Call:
## lm(formula = CM1$y ~ CM1$x)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.1653 -2.5817 -0.0381  2.2415  5.5486 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   7.6573     1.5988    4.79 0.000147 ***
## CM1$x         1.7725     0.1548   11.45 1.07e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.192 on 18 degrees of freedom
## Multiple R-squared:  0.8792, Adjusted R-squared:  0.8725 
## F-statistic:   131 on 1 and 18 DF,  p-value: 1.074e-09

names(antoine) #Les éléments du modèle 'antoine'.

##  [1] "coefficients"  "residuals"     "effects"       "rank"         
##  [5] "fitted.values" "assign"        "qr"            "df.residual"  
##  [9] "xlevels"       "call"          "terms"         "model"

*** Fin de la présentation***

Send an E-mail to Omran ALLATIF

CIRI - Centre International de Recherche en Infectiologie

Download this presentation