Les dépenses pendant le réveillon 2015
SKANJI-MOEZ
partie descriptive
importation de la base de donnée
bd<- read.delim("final.txt")
head(bd,10)## sexe Age Votre.occupation..
## 1 homme 20.0 Etudiant /élève.
## 2 homme 20.0 Etudiant /élève.
## 3 homme 32.5 Employé/ouvrier.
## 4 homme 32.5 Cadre supérieur /profession libérale.
## 5 homme 20.0 Etudiant /élève.
## 6 homme 32.5 Agriculteur/artisan/commerçant.
## 7 femme 20.0 Etudiant /élève.
## 8 homme 32.5 Agriculteur/artisan/commerçant.
## 9 femme 47.5 Sans emploi.
## 10 homme 20.0 Etudiant /élève.
## Vous.avez..passé.le.dernier.réveillon..
## 1 A l<U+0092>hôtel, restaurant, lounge , discothèque
## 2 A l<U+0092>hôtel, restaurant, lounge , discothèque
## 3 Chez vous
## 4 Chez vous
## 5 A l<U+0092>hôtel, restaurant, lounge , discothèque
## 6 Chez la famille ou les amis
## 7 A l<U+0092>étranger
## 8 Chez la famille ou les amis
## 9 Chez la famille ou les amis
## 10 A l<U+0092>étranger
## Quel.budget.en.dinars.avez.vous.consacré..pour..le.dernier.réveillon..
## 1 200
## 2 400
## 3 25
## 4 25
## 5 200
## 6 25
## 7 75
## 8 75
## 9 75
## 10 400
## Pendant.la.soirée.du..réveillon..vous.avez.consommé..
## 1 Uniquement des boissons sans alcool
## 2 Des boissons alcoolisées
## 3 Uniquement des boissons sans alcool
## 4 Uniquement des boissons sans alcool
## 5 Des boissons alcoolisées
## 6 Uniquement des boissons sans alcool
## 7 Uniquement des boissons sans alcool
## 8 Uniquement des boissons sans alcool
## 9 Uniquement des boissons sans alcool
## 10 Uniquement des boissons sans alcoolsummary(bd)## sexe Age Votre.occupation..
## femme:50 Min. :20.00 Agriculteur/artisan/commerçant. :10
## homme:50 1st Qu.:20.00 Cadre supérieur /profession libérale.:18
## Median :20.00 Employé/ouvrier. : 6
## Mean :33.17 Etudiant /élève. :49
## 3rd Qu.:47.50 Sans emploi. :17
## Max. :62.50
## Vous.avez..passé.le.dernier.réveillon..
## A l<U+0092>étranger : 8
## A l<U+0092>hôtel, restaurant, lounge , discothèque:11
## Chez la famille ou les amis :43
## Chez vous :38
##
##
## Quel.budget.en.dinars.avez.vous.consacré..pour..le.dernier.réveillon..
## Min. : 25.00
## 1st Qu.: 25.00
## Median : 75.00
## Mean : 99.25
## 3rd Qu.: 75.00
## Max. :400.00
## Pendant.la.soirée.du..réveillon..vous.avez.consommé..
## Des boissons alcoolisées :23
## Uniquement des boissons sans alcool:77
##
##
##
## colnames(bd)<-c("sexe","age","occupation","lieu", "budget","consommation")
table(bd$sexe)##
## femme homme
## 50 50df<-data.frame(group=c("Homme", "Femme") , value = c(50,50))
head(df)## group value
## 1 Homme 50
## 2 Femme 50pie chart : groupe de sexe
library(ggplot2)
library(scales)
bp<- ggplot(data=df, aes(x="", y=value, fill=group))+
geom_bar(width = 1, stat = "identity")
blank_theme <- theme_minimal()+
theme(
axis.title.x = element_blank(),
axis.title.y = element_blank(),
panel.border = element_blank(),
panel.grid=element_blank(),
axis.ticks = element_blank(),
plot.title=element_text(size=14, face="bold")
)
pie <- bp + coord_polar("y", start=0)
pie +scale_fill_brewer(palette="Blues")+theme_minimal()+blank_theme +
theme(axis.text.x=element_blank())+
geom_text(aes(y = value/3 + c(0, cumsum(value)[-length(value)]),
label = percent(value/100)), size=5)
groupement d’age
summary(bd$age)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 20.00 20.00 20.00 33.17 47.50 62.50table(bd$age)##
## 20 32.5 47.5 62.5
## 51 19 13 17df<-data.frame(group=c("[15,25[", "[25,40[","[40,55[","[55,70[") , value = c(51,19,13,17))
head(df)## group value
## 1 [15,25[ 51
## 2 [25,40[ 19
## 3 [40,55[ 13
## 4 [55,70[ 17bp<- ggplot(data=df, aes(x="", y=value, fill=group))+
geom_bar(width = 1, stat = "identity")
blank_theme <- theme_minimal()+
theme(
axis.title.x = element_blank(),
axis.title.y = element_blank(),
panel.border = element_blank(),
panel.grid=element_blank(),
axis.ticks = element_blank(),
plot.title=element_text(size=14, face="bold")
)
pie <- bp + coord_polar("y", start=0)
pie +scale_fill_brewer(palette="Blues")+theme_minimal()+blank_theme +
theme(axis.text.x=element_blank())+
geom_text(aes(y = value/3 + c(0, cumsum(value)[-length(value)]),
label = percent(value/100)), size=5) 
les occupations
table(bd$occupation)##
## Agriculteur/artisan/commerçant.
## 10
## Cadre supérieur /profession libérale.
## 18
## Employé/ouvrier.
## 6
## Etudiant /élève.
## 49
## Sans emploi.
## 17df<-data.frame(group=c("Etudiant /élève","Employé/ouvrier","Cadre supérieur
/profession libérale","Agriculteur/artisan/commerçant","Sans emploi") ,
value = c(49,6,18,10,17))
head(df)## group value
## 1 Etudiant /élève 49
## 2 Employé/ouvrier 6
## 3 Cadre supérieur \n /profession libérale 18
## 4 Agriculteur/artisan/commerçant 10
## 5 Sans emploi 17bp<- ggplot(data=df, aes(x="", y=value, fill=group))+
geom_bar(width = 1, stat = "identity")
bp
blank_theme <- theme_minimal()+
theme(
axis.title.x = element_blank(),
axis.title.y = element_blank(),
panel.border = element_blank(),
panel.grid=element_blank(),
axis.ticks = element_blank(),
plot.title=element_text(size=14, face="bold")
)
pie <- bp + coord_polar("y", start=0)
pie +scale_fill_brewer(palette="Dark2")+theme_minimal()+blank_theme +
theme(axis.text.x=element_blank())+
geom_text(aes(y = value/3 + c(0, cumsum(value)[-length(value)]),
label = percent(value/100)), size=5) 
les lieux de la soirée
table(bd$lieu)##
## A l<U+0092>étranger
## 8
## A l<U+0092>hôtel, restaurant, lounge , discothèque
## 11
## Chez la famille ou les amis
## 43
## Chez vous
## 38df<-data.frame(group=c("A l'hôtel, restaurant, lounge , discothèque"
,"Chez vous","A l'étranger","Chez la famille ou les amis")
, value = c(11,38,8,43))
head(df)## group value
## 1 A l'hôtel, restaurant, lounge , discothèque 11
## 2 Chez vous 38
## 3 A l'étranger 8
## 4 Chez la famille ou les amis 43bp<- ggplot(data=df, aes(x="", y=value, fill=group))+
geom_bar(width = 1, stat = "identity")
bp
blank_theme <- theme_minimal()+
theme(
axis.title.x = element_blank(),
axis.title.y = element_blank(),
panel.border = element_blank(),
panel.grid=element_blank(),
axis.ticks = element_blank(),
plot.title=element_text(size=14, face="bold")
)
pie <- bp + coord_polar("y", start=0)
pie +scale_fill_brewer(palette="Blues")+theme_minimal()+blank_theme +
theme(axis.text.x=element_blank())+
geom_text(aes(y = value/3 + c(0, cumsum(value)[-length(value)]),
label = percent(value/100)), size=5) 
la consommation
table(bd$consommation)##
## Des boissons alcoolisées Uniquement des boissons sans alcool
## 23 77df<-data.frame(group=c("des boissons alcolisées","des boissons vierges ") , value = c(23,77))
head(df)## group value
## 1 des boissons alcolisées 23
## 2 des boissons vierges 77bp<- ggplot(data=df, aes(x="", y=value, fill=group))+
geom_bar(width = 1, stat = "identity")
bp
blank_theme <- theme_minimal()+
theme(
axis.title.x = element_blank(),
axis.title.y = element_blank(),
panel.border = element_blank(),
panel.grid=element_blank(),
axis.ticks = element_blank(),
plot.title=element_text(size=14, face="bold")
)
pie <- bp + coord_polar("y", start=0)
pie +scale_fill_brewer(palette="Blues")+theme_minimal()+blank_theme +
theme(axis.text.x=element_blank())+
geom_text(aes(y = value/3 + c(0, cumsum(value)[-length(value)]),
label = percent(value/100)), size=5) 
#boxplot ##age&budget
p<-ggplot(bd, aes(factor(age),budget) )
p+geom_boxplot()
p<-ggplot(bd, aes(factor(sexe),budget) )+geom_boxplot(outlier.size = 3)+coord_flip()
p
#epartition de l’age selon le sexe
p <- ggplot(bd, aes(factor(sexe),age)) + geom_boxplot()+ geom_jitter()
p
opaa<-bd$age[which(bd$sexe == "homme")]
opaa## [1] 20.0 20.0 32.5 32.5 20.0 32.5 32.5 20.0 32.5 20.0 32.5 32.5 47.5 47.5
## [15] 32.5 20.0 20.0 20.0 20.0 20.0 20.0 20.0 20.0 32.5 20.0 20.0 32.5 32.5
## [29] 20.0 20.0 20.0 20.0 32.5 62.5 32.5 62.5 47.5 20.0 47.5 20.0 47.5 62.5
## [43] 62.5 62.5 62.5 62.5 62.5 47.5 62.5 62.5summary(opaa)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 20.00 20.00 32.50 35.05 47.50 62.50oppa2<-bd$age[which(bd$sexe == "femme")]
summary(oppa2)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 20.0 20.0 20.0 31.3 47.5 62.5boxplot du budget selon le sexe
p <- ggplot(bd, aes(factor(sexe),budget)) + geom_boxplot()+ geom_jitter()
p
opaa<-bd$budget[which(bd$sexe == "homme")]
opaa## [1] 200 400 25 25 200 25 75 400 200 75 25 75 75 75 25 25 75
## [18] 25 75 25 25 200 25 25 25 25 25 25 75 25 75 75 25 75
## [35] 200 400 75 75 75 25 25 25 25 25 75 75 75 200 200 200summary(opaa)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 25.0 25.0 75.0 92.5 75.0 400.0oppa2<-bd$budget[which(bd$sexe == "femme")]
summary(oppa2)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 25 25 75 106 75 400histogramme
qplot(factor(budget) , data=bd , geom="bar", fill=factor(sexe))
qplot(factor(budget) , data=bd , geom="bar", fill=factor(lieu))
qplot(factor(budget) , data=bd , geom="bar", fill=factor(occupation))
qplot(factor(budget) , data=bd , geom="bar", fill=factor(consommation))
sas
#intialisation
ypesr<-vector(length=20 , "numeric")
ypesr## [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0vi<-bd$budget
ech<-NULL
for (i in 1:20 ){
ech<-cbind(ech,sample(vi,20))
ypesr[i]<-mean(ech[,i])
}
ymoy<-mean(ypesr)
ymoy ## [1] 92.1875ic<-matrix(0,20,2)
sc<-vector(length=20,"numeric")
for(i in 1:20 ){
sc[i]<-(1/19)*sum((ech[,i]-ypesr[i])^2)
variance_pesr_estime<-0.04*sc
ic[i,]=c( ypesr[i]-1.96*sqrt(variance_pesr_estime[i]/20)
,ypesr[i]+1.96*sqrt(variance_pesr_estime[i]/20) )
langeur<-ic[,2]-ic[,1]
}
plot(c(1:20),ic[,1],ylim=c(min(min(ic[,2]),min(ic[,1]))
,max(max(ic[,2]),max(ic[,1]))),type="o",xaxt =
"n",xlab="les échantillons",ylab="Intervalle de confiance",col="blue")
par(new=TRUE)
plot(c(1:20),ic[,2],ylim=c(min(min(ic[,2]),min(ic[,1]))
,max(max(ic[,2]),max(ic[,1]))),type="o",xaxt =
"n",xlab="les échantillons",ylab="Intervalle de confiance",cex=1.5,font=3,col="red")
axis(1,at=c(1:20))
moy=mean(ypesr[i])
TAB<-cbind(ypesr,variance_pesr_estime, ic[,1],ic[,2],langeur)
TAB## ypesr variance_pesr_estime langeur
## [1,] 123.75 540.72368 113.55872 133.94128 20.382562
## [2,] 100.00 526.31579 89.94541 110.05459 20.109176
## [3,] 112.50 559.21053 102.13597 122.86403 20.728064
## [4,] 78.75 340.19737 70.66637 86.83363 16.167265
## [5,] 75.00 184.21053 69.05163 80.94837 11.896749
## [6,] 118.75 570.72368 108.27982 129.22018 20.940354
## [7,] 96.25 377.03947 87.73990 104.76010 17.020193
## [8,] 92.50 345.00000 84.35951 100.64049 16.280983
## [9,] 97.50 536.57895 87.34785 107.65215 20.304294
## [10,] 102.50 397.10526 93.76639 111.23361 17.467224
## [11,] 93.75 504.93421 83.90176 103.59824 19.696473
## [12,] 72.50 189.21053 66.47144 78.52856 12.057124
## [13,] 62.50 322.36842 54.63104 70.36896 15.737919
## [14,] 102.50 515.52632 92.54901 112.45099 19.901989
## [15,] 97.50 536.57895 87.34785 107.65215 20.304294
## [16,] 72.50 307.63158 64.81301 80.18699 15.373988
## [17,] 105.00 385.78947 96.39172 113.60828 17.216555
## [18,] 53.75 71.77632 50.03694 57.46306 7.426115
## [19,] 126.25 853.88158 113.44322 139.05678 25.613557
## [20,] 60.00 114.21053 55.31625 64.68375 9.367509colnames(TAB)<-c("Y_PESR" , "var_PESR", "IC_PESR_MIN" , "IC_PESR_MAX" ,"longueur_IC")
TAB## Y_PESR var_PESR IC_PESR_MIN IC_PESR_MAX longueur_IC
## [1,] 123.75 540.72368 113.55872 133.94128 20.382562
## [2,] 100.00 526.31579 89.94541 110.05459 20.109176
## [3,] 112.50 559.21053 102.13597 122.86403 20.728064
## [4,] 78.75 340.19737 70.66637 86.83363 16.167265
## [5,] 75.00 184.21053 69.05163 80.94837 11.896749
## [6,] 118.75 570.72368 108.27982 129.22018 20.940354
## [7,] 96.25 377.03947 87.73990 104.76010 17.020193
## [8,] 92.50 345.00000 84.35951 100.64049 16.280983
## [9,] 97.50 536.57895 87.34785 107.65215 20.304294
## [10,] 102.50 397.10526 93.76639 111.23361 17.467224
## [11,] 93.75 504.93421 83.90176 103.59824 19.696473
## [12,] 72.50 189.21053 66.47144 78.52856 12.057124
## [13,] 62.50 322.36842 54.63104 70.36896 15.737919
## [14,] 102.50 515.52632 92.54901 112.45099 19.901989
## [15,] 97.50 536.57895 87.34785 107.65215 20.304294
## [16,] 72.50 307.63158 64.81301 80.18699 15.373988
## [17,] 105.00 385.78947 96.39172 113.60828 17.216555
## [18,] 53.75 71.77632 50.03694 57.46306 7.426115
## [19,] 126.25 853.88158 113.44322 139.05678 25.613557
## [20,] 60.00 114.21053 55.31625 64.68375 9.367509avec remise
set.seed(100)
ic=matrix(0,20,2)
vari=vector(length=20,"numeric")
ypear=vector(length=20,"numeric")
ech=NULL
for (i in 1:20){
ech<-cbind(ech,sample(vi,20,replace=TRUE) )
ypear[i]=mean(ech[,i])
}
ymoy=mean(ypear)
ymoy## [1] 97.0625s=vector("length"=20,"numeric")
for(i in 1:20)
{s[i]=(1/19)*sum((ech[,i]-ypear[i])^2)}
variance_pear_estim<-s/20
variance_pear_estim## [1] 657.89474 703.61842 86.43092 1172.36842 461.84211 764.72039
## [7] 896.29934 721.95724 896.29934 191.69408 665.37829 642.35197
## [13] 458.79934 271.71053 425.24671 665.37829 532.48355 470.72368
## [19] 665.37829 147.69737for(i in 1:20){
ic[i,]=c(ypear[i]-1.96*sqrt(variance_pear_estim[i]/20),ypear[i]
+1.96*sqrt(variance_pear_estim[i]/20))}
ic## [,1] [,2]
## [1,] 88.75863 111.24137
## [2,] 115.87455 139.12545
## [3,] 54.67549 62.82451
## [4,] 144.99372 175.00628
## [5,] 75.58137 94.41863
## [6,] 99.13028 123.36972
## [7,] 100.62898 126.87102
## [8,] 86.97403 110.52597
## [9,] 100.62898 126.87102
## [10,] 60.18200 72.31800
## [11,] 97.44487 120.05513
## [12,] 80.14221 102.35779
## [13,] 89.36245 108.13755
## [14,] 87.77572 102.22428
## [15,] 69.71222 87.78778
## [16,] 97.44487 120.05513
## [17,] 71.13667 91.36333
## [18,] 72.99124 92.00876
## [19,] 97.44487 120.05513
## [20,] 47.17368 57.82632longeur_pear<-ic[,2]-ic[,1]
longeur_pear## [1] 22.482742 23.250895 8.149025 30.012566 18.837264 24.239430 26.242041
## [8] 23.551947 26.242041 12.135996 22.610251 22.215577 18.775109 14.448551
## [15] 18.075551 22.610251 20.226660 19.017529 22.610251 10.652645plot(c(1:20),ic[,1],ylim=c(min(min(ic[,2]),min(ic[,1])),max(max(ic[,2]),max(ic[,1]))),
type="o",xaxt
= "n",xlab="les échantillons",ylab="Intervalle de confiance",col="blue")
par(new=TRUE)
plot(c(1:20),ic[,2],ylim=c(min(min(ic[,2]),min(ic[,1])),max(max(ic[,2])
,max(ic[,1]))),type="o",xaxt
= "n",xlab="les échantillons",ylab="Intervalle de confiance",cex=1.5,font=3,col="red")
axis(1,at=c(1:20))
TAB2<-cbind(ypear,variance_pear_estim , ic [,1] , ic[,2] , longeur_pear)
TAB2## ypear variance_pear_estim longeur_pear
## [1,] 100.00 657.89474 88.75863 111.24137 22.482742
## [2,] 127.50 703.61842 115.87455 139.12545 23.250895
## [3,] 58.75 86.43092 54.67549 62.82451 8.149025
## [4,] 160.00 1172.36842 144.99372 175.00628 30.012566
## [5,] 85.00 461.84211 75.58137 94.41863 18.837264
## [6,] 111.25 764.72039 99.13028 123.36972 24.239430
## [7,] 113.75 896.29934 100.62898 126.87102 26.242041
## [8,] 98.75 721.95724 86.97403 110.52597 23.551947
## [9,] 113.75 896.29934 100.62898 126.87102 26.242041
## [10,] 66.25 191.69408 60.18200 72.31800 12.135996
## [11,] 108.75 665.37829 97.44487 120.05513 22.610251
## [12,] 91.25 642.35197 80.14221 102.35779 22.215577
## [13,] 98.75 458.79934 89.36245 108.13755 18.775109
## [14,] 95.00 271.71053 87.77572 102.22428 14.448551
## [15,] 78.75 425.24671 69.71222 87.78778 18.075551
## [16,] 108.75 665.37829 97.44487 120.05513 22.610251
## [17,] 81.25 532.48355 71.13667 91.36333 20.226660
## [18,] 82.50 470.72368 72.99124 92.00876 19.017529
## [19,] 108.75 665.37829 97.44487 120.05513 22.610251
## [20,] 52.50 147.69737 47.17368 57.82632 10.652645colnames(TAB2)<-cbind("Y_PEAR" , "var_PEAR" , "IC_PEAR_MIN" , "IC_PEAR_MAX" , "longueur_pear")
TAB2## Y_PEAR var_PEAR IC_PEAR_MIN IC_PEAR_MAX longueur_pear
## [1,] 100.00 657.89474 88.75863 111.24137 22.482742
## [2,] 127.50 703.61842 115.87455 139.12545 23.250895
## [3,] 58.75 86.43092 54.67549 62.82451 8.149025
## [4,] 160.00 1172.36842 144.99372 175.00628 30.012566
## [5,] 85.00 461.84211 75.58137 94.41863 18.837264
## [6,] 111.25 764.72039 99.13028 123.36972 24.239430
## [7,] 113.75 896.29934 100.62898 126.87102 26.242041
## [8,] 98.75 721.95724 86.97403 110.52597 23.551947
## [9,] 113.75 896.29934 100.62898 126.87102 26.242041
## [10,] 66.25 191.69408 60.18200 72.31800 12.135996
## [11,] 108.75 665.37829 97.44487 120.05513 22.610251
## [12,] 91.25 642.35197 80.14221 102.35779 22.215577
## [13,] 98.75 458.79934 89.36245 108.13755 18.775109
## [14,] 95.00 271.71053 87.77572 102.22428 14.448551
## [15,] 78.75 425.24671 69.71222 87.78778 18.075551
## [16,] 108.75 665.37829 97.44487 120.05513 22.610251
## [17,] 81.25 532.48355 71.13667 91.36333 20.226660
## [18,] 82.50 470.72368 72.99124 92.00876 19.017529
## [19,] 108.75 665.37829 97.44487 120.05513 22.610251
## [20,] 52.50 147.69737 47.17368 57.82632 10.652645PISR
library(sampling)
y=bd[,5]## Variable d'interet
x=bd$age## Variable auxiliaire
x=na.omit(x)
n=20
x## [1] 20.0 20.0 32.5 32.5 20.0 32.5 20.0 32.5 47.5 20.0 20.0 20.0 20.0 20.0
## [15] 20.0 20.0 20.0 32.5 20.0 20.0 20.0 20.0 20.0 20.0 20.0 20.0 32.5 32.5
## [29] 47.5 32.5 47.5 32.5 20.0 20.0 20.0 20.0 20.0 20.0 20.0 32.5 20.0 20.0
## [43] 47.5 20.0 20.0 32.5 20.0 20.0 20.0 20.0 32.5 20.0 20.0 20.0 20.0 20.0
## [57] 20.0 32.5 32.5 32.5 20.0 20.0 32.5 20.0 20.0 20.0 32.5 62.5 32.5 47.5
## [71] 62.5 47.5 20.0 20.0 62.5 47.5 47.5 20.0 32.5 47.5 20.0 47.5 62.5 62.5
## [85] 62.5 62.5 62.5 62.5 62.5 62.5 62.5 62.5 62.5 62.5 47.5 47.5 47.5 62.5
## [99] 62.5 20.0echantillon=matrix(0,nrow=20,ncol=20)## initialisation de la matrice contenant 20*20
IC=matrix(0,nrow=20,ncol=2) ## initialisation de la matrice contenant 20 IC
ybarPISR=NULL ## initialisation du vecteur contenant les 20 estimateurs de
varPISR=NULL ## initialisation du vecteur contenant 20 estimateurs de la variance
pik=inclusionprobabilities(x,n) ## les probabilites d'inclusion
pik## [1] 0.1205727 0.1205727 0.1959307 0.1959307 0.1205727 0.1959307 0.1205727
## [8] 0.1959307 0.2863602 0.1205727 0.1205727 0.1205727 0.1205727 0.1205727
## [15] 0.1205727 0.1205727 0.1205727 0.1959307 0.1205727 0.1205727 0.1205727
## [22] 0.1205727 0.1205727 0.1205727 0.1205727 0.1205727 0.1959307 0.1959307
## [29] 0.2863602 0.1959307 0.2863602 0.1959307 0.1205727 0.1205727 0.1205727
## [36] 0.1205727 0.1205727 0.1205727 0.1205727 0.1959307 0.1205727 0.1205727
## [43] 0.2863602 0.1205727 0.1205727 0.1959307 0.1205727 0.1205727 0.1205727
## [50] 0.1205727 0.1959307 0.1205727 0.1205727 0.1205727 0.1205727 0.1205727
## [57] 0.1205727 0.1959307 0.1959307 0.1959307 0.1205727 0.1205727 0.1959307
## [64] 0.1205727 0.1205727 0.1205727 0.1959307 0.3767898 0.1959307 0.2863602
## [71] 0.3767898 0.2863602 0.1205727 0.1205727 0.3767898 0.2863602 0.2863602
## [78] 0.1205727 0.1959307 0.2863602 0.1205727 0.2863602 0.3767898 0.3767898
## [85] 0.3767898 0.3767898 0.3767898 0.3767898 0.3767898 0.3767898 0.3767898
## [92] 0.3767898 0.3767898 0.3767898 0.2863602 0.2863602 0.2863602 0.3767898
## [99] 0.3767898 0.1205727#tirage des echantillons:
for (i in 1:20){
s=UPsystematic(pik)
echantillon[,i]=y[s==1][1:20]
ybarPISR[i]=mean(echantillon[,i]) ##les estimateurs de la moyenne
varPISR[i]=var(echantillon[,i])} ##les estimateurs de la variance
moy_PISR=mean(ybarPISR)
#calcul des intervalles de confiance:
for(i in 1:20){
IC[i,]=c(ybarPISR[i]-1.96*sqrt(var(ybarPISR)/20), ##borne inf IC
ybarPISR[i]+1.96*sqrt(var(ybarPISR)/20))} ##borne sup IC
longueurIC=IC[,2]-IC[,1]
#representation des IC:
plot(c(1:20),IC[,1],col="deepskyblue",ylim=c(min(min(IC[,2]),
min(IC[,1])),
max(max(IC[,2]),max(IC[,1]))))
icPISRmin=IC[,1]
icPISRmax=IC[,2]
lines(seq(1,20,by=1),IC[,1],col="deepskyblue")
points(seq(1,20,by=1),IC[,2],col="deeppink")
lines(seq(1,20,by=1),IC[,2],col="deeppink")
TAB3=cbind(ybarPISR,varPISR,icPISRmin,icPISRmax,longueurIC,moy_PISR)
TAB3## ybarPISR varPISR icPISRmin icPISRmax longueurIC moy_PISR
## [1,] 102.50 9927.632 91.63524 113.36476 21.72953 103.8125
## [2,] 116.25 12320.724 105.38524 127.11476 21.72953 103.8125
## [3,] 117.50 13296.053 106.63524 128.36476 21.72953 103.8125
## [4,] 121.25 13899.671 110.38524 132.11476 21.72953 103.8125
## [5,] 68.75 3741.776 57.88524 79.61476 21.72953 103.8125
## [6,] 100.00 10197.368 89.13524 110.86476 21.72953 103.8125
## [7,] 77.50 4467.105 66.63524 88.36476 21.72953 103.8125
## [8,] 77.50 9730.263 66.63524 88.36476 21.72953 103.8125
## [9,] 77.50 9730.263 66.63524 88.36476 21.72953 103.8125
## [10,] 117.50 13296.053 106.63524 128.36476 21.72953 103.8125
## [11,] 143.75 20518.092 132.88524 154.61476 21.72953 103.8125
## [12,] 96.25 9425.987 85.38524 107.11476 21.72953 103.8125
## [13,] 125.00 12171.053 114.13524 135.86476 21.72953 103.8125
## [14,] 117.50 18559.211 106.63524 128.36476 21.72953 103.8125
## [15,] 143.75 20518.092 132.88524 154.61476 21.72953 103.8125
## [16,] 120.00 18197.368 109.13524 130.86476 21.72953 103.8125
## [17,] 128.75 17978.618 117.88524 139.61476 21.72953 103.8125
## [18,] 77.50 4467.105 66.63524 88.36476 21.72953 103.8125
## [19,] 85.00 9236.842 74.13524 95.86476 21.72953 103.8125
## [20,] 62.50 2796.053 51.63524 73.36476 21.72953 103.8125PIAR
y=bd$budget
n=20
pk=inclusionprobabilities(y,n)/n
sum(pk)## [1] 1v=c(0,cumsum(pk))
ech=matrix(0,nrow=20,ncol=20)
IC=matrix(0,nrow=20,ncol=2)
y.bar.piar=NULL
var.piar=NULL
##tirage
##Tirage de 20 échantillons
for (i in 1:20){
s=NULL
for(j in 1:20){
u=runif(1,0,1)
for(k in 2:101){
if((v[k-1]<=u)&(v[k]>u)) {
s=c(s,k-1)
} }}
ech[,i]=y[s]
y.bar.piar[i]=mean(ech[,i])
var.piar[i]=var(ech[,i])
}
moy_y.bar.piar=mean(y.bar.piar)
##Intervalles de confiance
for(i in 1:20)
{IC[i,]=c(y.bar.piar[i]-1.96*sqrt(var(y.bar.piar)/20),y.bar.piar[i]+1.96*sqrt(var(y.bar.piar)/20))
}
longueur=IC[,2]-IC[,1]
moy_y.bar.piar=mean(y.bar.piar)
##Représentation des IC
plot(c(1:20),IC[,1],ylim=c(min(min(IC[,2]),min(IC[,1])),max(max(IC[,2]),max(IC[,1]))),type="o",xaxt =
"n",xlab="les échantillons",ylab="Intervalle de confiance",col="blue")
par(new=TRUE)
plot(c(1:20),IC[,2],ylim=c(min(min(IC[,2]),min(IC[,1])),max(max(IC[,2]),max(IC[,1]))),type="o",xaxt =
"n",xlab="les échantillons",ylab="Intervalle de confiance",cex=1.5,font=3,col="red")
axis(1,at=c(1:20))
piar=cbind(y.bar.piar,var.piar,IC[,2],IC[,1] , longueur)
colnames(piar)<-c("y.bar.piar" , "var.piar" , "icPIARmax" , "icPIARmin" , "longueur_PIAR")
piar## y.bar.piar var.piar icPIARmax icPIARmin longueur_PIAR
## [1,] 230.00 19513.16 245.1726 214.8274 30.34525
## [2,] 176.25 17597.04 191.4226 161.0774 30.34525
## [3,] 276.25 21939.14 291.4226 261.0774 30.34525
## [4,] 281.25 16570.72 296.4226 266.0774 30.34525
## [5,] 231.25 23149.67 246.4226 216.0774 30.34525
## [6,] 212.50 23519.74 227.6726 197.3274 30.34525
## [7,] 223.75 17728.62 238.9226 208.5774 30.34525
## [8,] 246.25 19491.78 261.4226 231.0774 30.34525
## [9,] 198.75 18649.67 213.9226 183.5774 30.34525
## [10,] 251.25 25623.36 266.4226 236.0774 30.34525
## [11,] 212.50 23519.74 227.6726 197.3274 30.34525
## [12,] 276.25 21939.14 291.4226 261.0774 30.34525
## [13,] 266.25 21333.88 281.4226 251.0774 30.34525
## [14,] 210.00 16144.74 225.1726 194.8274 30.34525
## [15,] 331.25 12952.30 346.4226 316.0774 30.34525
## [16,] 226.25 25031.25 241.4226 211.0774 30.34525
## [17,] 237.50 21875.00 252.6726 222.3274 30.34525
## [18,] 222.50 25256.58 237.6726 207.3274 30.34525
## [19,] 251.25 28583.88 266.4226 236.0774 30.34525
## [20,] 236.25 18254.93 251.4226 221.0774 30.34525