Tirage aléatoire dans une matrice et simulations

Version imprimable

Bonjour

Je me suis inscrit il y a peu et c’est donc ma première question. J’espère que tout sera clair, si non, je ferais en sorte d’apporter les informations manquantes.

Je cherche à réaliser un test de puissance afin de savoir le nombre d’échantillons nécessaires pour mon échantillonnage (Intervalle de confiance correct…)

1ère partie :

• Je possède un jeu de données sous forme de matrice de 20 lignes et 15 colonnes (noté a=read.table…)

http://img11.hostingpics.net/pics/724722tab.jpg

• Je cherche à réaliser un tirage aléatoire de X valeurs dans des sous-unités (5*5, soit 12 sous-unités). Ne sachant pas comment faire pour que, dans la matrice principale, R comprennent qu’il y a 12 sous-unités, j’ai créé les sous-unités de la manière suivante :
Code:

1 2 #Définition des sous matrices A1=a[1:5,1:5];A2=a[6:10,1:5];A3=a[11:15,1:5];A4=a[16:20,1:5];A5=a[1:5,6:10];A6=a[6:10,6:10];A7=a[11:15,6:10];A8=a[16:20,6:10];A9=a[1:5,11:15];A10=a[6:10,11:15];A11=a[11:15,11:15];A12=a[16:20,11:15]
• A priori, je ne peux pas utiliser la fonction sample comme ça, je dois
Code:

1 2 3 4 transformer ces "sous matrices" en réelles matrices : A1 <- as.matrix(A1);A2 <- as.matrix(A2);A3 <- as.matrix(A3);A4 <- as.matrix(A4); A5 <- as.matrix(A5);A6 <- as.matrix(A6);A7 <- as.matrix(A7);A8 <- as.matrix(A8); A9 <- as.matrix(A9);A10 <- as.matrix(A10);A11 <- as.matrix(A11);A12 <- as.matrix(A12);
et là je peux réaliser la fonction sample, c'est bien cela ? (je bloquais avant parce que je n'utilisais pas as.matrix et il ne devait pas comprendre que c'était une matrice vu qu'il me sélectionnait toute la colonne et non une seule valeur).

2ème partie :
• Là, la fonction sample fonctionne. Je l'utilise en définissant

Code:

samp<-list(s1=sample(A1,3,replace=F),sample(A2,3,replace=F),...,s12=sample(A12,3,replace=F))

J'ai donc une liste avec mes 12 strates composées chacune d'une valeur.
Je voudrais sélectionner uniquement celles où il y a au moins une valeur positive. J'avais pensé à utiliser la fonction if, mais je ne sais pas trop comment lui dire "si au moins une valeur est positive, je garde le lot, sinon je ne le prends pas en compte"

Ceci me permettrait par la suite de réaliser la formule suivante : (pour exemple : 3 strates sélectionnées)
Code:

1 2 3 4 5 6 7 8 9 10 11 12 13 N <- 12 #Nombre de strates Nh <-25 #Nombre de mailles par strate nh <- X #Nombre de mailles secondaires par strate v <- 3 #nombre de strates où l’espèce est observée lors de la première étape d’échantillonnage M <- v * Nh + (N-v) * nh #nombre de mailles échantillonnées mh <- c(7,14,14) #nombre de mailles dans le réseau (non vides dans la strate) y <- c(15,66,73) #valeurs de la variable d’intérêt observée pik <- 1 - choose(Nh-mh,3)/choose(Nh,nh) tauhat <- sum(y/pik) #estimation de la population totale (HT) var.tauhat <- sum(((1-pik)*y^2)/(pik^2)) c(tauhat,sqrt(var.tauhat)) #Estimation du total et Ecart-type IC95<- c(tauhat-((1.96*sqrt(var.tauhat))/sqrt(M)), tauhat+((1.96*sqrt(var.tauhat))/sqrt(M))) IC99<- c(tauhat-((2.58 *sqrt(var.tauhat))/sqrt(M)), tauhat+((2.58 *sqrt(var.tauhat))/sqrt(M)))
Dois-je définir mh et y manuellement à chaque fois ou est-ce possible de lui indiquer que si l'unité A4 est sélectionnée (avec 1 valeur positive parmi les X sélectionnées) alors il y a mh=7 et y=66 ? (et donc définir chaque mh et y pour chaque strate)
Pareil pour v, puis-je lui indiquer que si une strate est sélectionnée, il me compte comme 1 et donc que v soit la somme des 1 ?

A terme je voudrais lancer la simulation un certain nombre de fois pour avoir différents intervalles de confiance pour différents nombres de sous-unités tirées aléatoirement (si je prends 3 valeurs dans la sous-unité ou 6 valeurs pour 100 simulations, je n’aurais pas le même intervalle de confiance et donc quel est le nombre adéquat à réaliser). Autrement dit faire 100 tirages aléatoires et calculer les IC pour chacun des tirages.

Voilà, j’espère avoir été assez clair, c’est un peu flou pour moi au final, surtout si je dois réaliser une boucle pour la simulation.

Je vous prie de m’excuser pour la non clarté si c’est le cas.

Merci d'avance,
Cordialement.

06/04/2015, 21h46
chautho

Bonjour,

Je vais déjà tâcher de répondre à un bout de ta question :

Citation:

Je voudrais sélectionner uniquement celles où il y a au moins une valeur positive. J'avais pensé à utiliser la fonction if, mais je ne sais pas trop comment lui dire "si au moins une valeur est positive, je garde le lot, sinon je ne le prends pas en compte"

Pour cela je prends ta liste "samp". Pour savoir si il y a au moins un nombre positif dans chaque élément de la liste ou non, je te conseille d'utiliser la fonction lapply().
Exemple :

Code:

lapply(samp,function(x) (sum(x > 0)) )

Ceci te renvoie un vecteur de la même longueur que ta liste en mettant un TRUE si il y a au moins un nombre positif, et un FALSE sinon.

En espérant t'avoir aidé, bon courage !
07/04/2015, 11h32
Geofv

Bonjour !

Merci pour la réponse.
Cependant, ça ne me donne que de l'information et non pas une sélection. Au lieu d'avoir une liste avec 3 chiffres (et les valeurs associées), j'obtiens une liste avec 1 chiffre (0 à 3) mais je souhaiterais qu'il ne me garde pas les listes où il y a la valeur 0.

En résumé :
j'ai [[1]] 0 (ex 0 0 0) [[2]] 2 (ex 0 5 1) ...
et je souhaiterais que ma liste ne contienne que les 1/2/3 et non pas les 0.

Au final, il ne me resterait que les [[X]] possédant au moins 1 valeur positive (disons 2 sur 12 par exemple).
Le but est de pouvoir l'implémenter dans une sorte de boucle afin de lancer la simulation plusieurs fois et donc ne pas avoir à sélectionner manuellement les listes ayant une valeur positive.
07/04/2015, 11h39
chautho

Bonjour,

Désolé mais je ne comprends vraiment pas ta question, aurais-tu l'exemple (en code R) d'une (ou plusieurs) liste(s) et ce que tu veux obtenir avec ?

Je ne comprends pas trop si tu veux retenir uniquement les listes ou supprimer les valeurs au sein des éléments de la liste...

Merci !
07/04/2015, 11h53
Geofv
Bonjour,

j'ai :

Code:

list(c(0,0,0),c(0,0,0),c(0,0,0),c(0,0,1),c(0,5,0),c(1,2,0),c(0,0,1),c(0,0,0),c(0,0,0),c(0,0,0),c(0,0,0),c(0,0,0))

Je voudrais avoir au final :
Code:

1 2 3 4 5 6 7 8 9 10 11 12 [[4]] [1] 0 0 1 [[5]] [1] 0 5 0 [[6]] [1] 1 2 0 [[7]] [1] 0 0 1
Comme cela, je n'aurais que les "listes" qui m'intéressent et donc je pourrais associer la valeur qui y est associée (par exemple pour [[4]] :mh4=1 ; y4=1 ; [[5]] mh5=7 ; y5=15 etc.) et l'idéal serait qu'en plus de me sélectionner seulement ces valeurs, une fois qu'elles sont sélectionnées, il m'y associe les deux valeurs correspondantes (mh et y).

Bonjour,

Avec ce que je t'ai donné ceci est relativement simple d'extraire les éléments d'une liste que tu veux !
Code:

1 2 3 4 5 6 liste <- list(c(0,0,0),c(0,0,0),c(0,0,0), c(0,0,1),c(0,5,0),c(1,2,0), c(0,0,1),c(0,0,0),c(0,0,0), c(0,0,0),c(0,0,0),c(0,0,0)) liste <- liste[as.logical(lapply(liste,function(x) sum(x>0) > 0 ))]
Voilà !

07/04/2015, 14h13
Geofv

Bonjour,

Merci, en effet ça fonctionne.
Il me reste donc à trouver comment lui dire que si [[4]] est sélectionné (valeur positive) alors mh4=X et y=Z (et faire ça de manière générale pour pouvoir le relancer plusieurs fois --> si c'est [[1]] qui a des valeurs positive alors mh1=.. et y1=.. , si c'est 2 alors .... etc)

merci
07/04/2015, 20h37
chautho

Essaye de ton coté, envoie nous tes codes si tu as des erreurs ou si tu es bloqué sur certain point :)

N'hésite par à mettre ton poste en résolu et nous donner la solution finale quand tu auras fini :)

1 pièce(s) jointe(s)

Bonjour,

Votre code m'a bien débloqué en fait, je pense avoir réussi ce que je voulais (bon peut être que c'est de la bidouille mais bon)
tout d'abord j'ai défini mes valeurs pour chaque mh et y (possible en faisant une somme de chaque pixel à priori, mais là je testais manuellement)

J'applique mon tirage aléatoire (n1=3)

Code:

samp2<-list(s1=sample(A1,n1,replace=F),s2=sample(A2,n1,replace=F),...,s12=sample(A12,n1,replace=F))

J'applique :

Code:

st <- lapply(samp2,function(x) (sum(x > 0)) )

pour savoir s'il y a des valeurs positives, et j'affecte mes valeurs mh aux valeurs

Code:

d2<-c(ifelse(st$s1>0,mh1,0),ifelse(st$s2>0 ,mh2,0),...,ifelse(st$s12>0 ,mh12,0)

de même pour les valeurs y

Code:

d3<- c(ifelse(d2[1]>0 ,0,0),ifelse(d2[2]>0 ,0,0),ifelse(d2[3]>0 ,0,0),...,ifelse(d2[12]>0 ,0,0)

Je ne sélectionne que les valeurs positives:
Code:

1 2 d2T <- d2[as.logical(lapply(d2,function(x) sum(x>0) > 0 ))] d3T <- d3[as.logical(lapply(d3,function(x) sum(x>0) > 0 ))]
Je compte le nombre d'unités qui présentent une valeur positive (au moins une parmi les 25 sous unités de la grande unité)
Code:

1 2 g<-c(1:12) v<-sum(ifelse(st[g]>0,1,0))
et après pour coller à ma formule:
Code:

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 mh<-c(d2T[g]) mh<-na.omit(mh) mh<-c(mh) y<-c(d3T[g]) y<-na.omit(y) y<-c(y) N <- 12 #Nombre de GU Nh <-25 #Nombre de mailles par strate nh <- n1 #Nombre de mailles secondaires par strate M <- v * Nh + (N-v) * nh #nombre de mailles échantillonnées pik <- 1 - choose(Nh-mh,3)/choose(Nh,nh) tauhat <- sum(y/pik) #estimation de la population totale (HT) var.tauhat <- sum(((1-pik)*y^2)/(pik^2)) c(tauhat,sqrt(var.tauhat)) # Estimation du total et Ecart type IC95<- c(tauhat-((1.96*sqrt(var.tauhat))/sqrt(M)), tauhat+((1.96*sqrt(var.tauhat))/sqrt(M)))
Maintenant, il faudrait que je réalise une boucle (for in 1: ..) pour pouvoir le réaliser plusieurs fois et avoir plusieurs valeurs et plusieurs intervalles de confiance et à terme avoir un graphique du genre : Pièce jointe 174072

et donc à partir de n1=3 j'ai un intervalle de confiance correct et donc je peux choisir 3 carrés comme valeur

Bonjour !

J'apporte un complément à ma question !
Finalement, une personne est venue m'apporter une vision différente et il a fallu que je retravaille sur l'ensemble du jeu de données.
En effet, je faisais les choses manuellement et donc pour pouvoir l'appliquer à un grand jeu de données ce n'était pas terrible.

Étant donné qu'il n'y a pas d'information de localisation du pixel, le groupe de 25 pixels à été mise dans une seule colonne et ainsi de suite pour chaque unité

ex :

Code:

Soit 300 lignes et 2 colonnées (col 1 = nom des unités, col2= effectif par pixel)

par la suite j'ai appliqué le code :

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
Tab=read.table("Test2.csv",sep=";",h=F)
 
#Somme des valeurs de chaque pixel pour une Grande Unité (GU)
#Donne les effectifs par GU (de 1 a 12)
library(doBy) #Appel de la librairie pour réaliser le summaryBy
Someff<-summaryBy(V2~V1,data=Tab,FUN=sum)
 
#changement des noms des colonnes de Someff
names(Someff)<- c("V1","V2.eff")
 
#Effectif total (toutes GU)
Efftot <- colSums(Tab["V2"])
 
#Nombre de cases non vides par GU
Pixnv<-data.frame(Someff["V1"],"V2"=(25-table(Tab)[,1]))
names(Pixnv)<- c("V1","V2.NV")
 
#Séparation de chaque Grande unité (liste)
Tab.split <- split(Tab,as.factor(Tab$V1))
 
#Nombre d'élements à tirer
n1=3
 
#Tirage aléatoire parmi chaque GU
Tab.sample <- lapply(Tab.split,function(x) x<-x[sample(1:nrow(x),n1,replace=F),])
 
#Transformation en dataframe des listes obtenues après tirage aléatoire
Tab.sample.df<-do.call(rbind.data.frame, Tab.sample) #Obtention data frame avec les 3 valeurs tirées par GU
 
#Affichage du nombre de valeur supérieures à "0" par GU
Tabs.df<-summaryBy(V2~V1,data=Tab.sample.df,FUN=function(x) (sum(x > 0)))
 
#Changement des noms des colonnes du data frame
names(Tabs.df)<- c("V1.1","V2.1")
attach(Tabs.df)
 
#Compilation des data frame m et Som (via noms des GU)
#nouveau data frame avec somme des effectifs par GU et nombre valeurs 
#positives avec tirage aléatoire de n1 valeurs
T1<-merge(Tabs.df, Someff , by.x = "V1.1", by.y = "V1")
 
#Compilation du data frame créé et data frame nombre pixels non vides par GU
#GU // Somme eff // Nombre de valeurs >0 parmi n1 // Nombre pixels non vides
Tab.merge<- merge(T1,Pixnv,by.x = "V1.1", by.y = "V1")
 
#Création data frame avec valeur de y (eff par GU)
d2<-data.frame(Someff["V1"],ifelse(Tab.merge[,2]>0, Tab.merge[,3], 0))
names(d2)<- c("V1","y")
 
#Création data frame avec valeur de mh (eff par GU)
d3<-data.frame(Someff["V1"],ifelse(Tab.merge[,2]>0, Tab.merge[,4], 0))
names(d3)<- c("V1","mh")
 
#Transformation des data frame en vecteur
d2.1<- as.vector(d2[,"y"])
d3.1<- as.vector(d3[,"mh"])
 
#Sélection des valeurs supérieures à 0 
y <- d2.1[as.logical(lapply(d2.1,function(x) sum(x>0) > 0 ))]
mh <- d3.1[as.logical(lapply(d3.1,function(x) sum(x>0) > 0 ))]
 
#Formule de la Variable v 
v<-sum(ifelse(Tab.merge["V2.1"]>0,1,0))
 
#Formule générale
N <- 12 #Nombre de GU 
Nh <-25 #Nombre de mailles par strate
nh <- n1 #Nombre de mailles secondaires par strate 
M <- v * Nh + (N-v) * nh #nombre de mailles échantillonnées
pik <- 1 - choose(Nh-mh,3)/choose(Nh,nh) 
tauhat <- sum(y/pik) #estimation de la population totale (HT)
var.tauhat <- sum(((1-pik)*y^2)/(pik^2))
c(tauhat,sqrt(var.tauhat)) # Estimation du total et Ecart type
IC95<- c(tauhat-((1.96*sqrt(var.tauhat))/sqrt(M)), tauhat+((1.96*sqrt(var.tauhat))/sqrt(M)))

Voilà, en espérant avoir pu aider quelqu'un également ;)

Geoffroy

Ps: cet exemple est réalisé dans le cadre d'un test pour un protocole d'échantillonnage nommé "Complete allocation sampling". Si vous souhaitez plus de détail sur ce protocole, n'hésitez pas, je vous ferai suivre quelques publications intéressantes.