Randomiser des conceptions expérimentales équilibrées

Question

J'écris du code pour générer des conceptions expérimentales équilibrées pour des études de marché, spécifiquement pour une utilisation dans l'analyse conjointe et l'échelle de différence maximale.

La première étape consiste à générer une conception de bloc incomplet partiellement équilibré (PBIB). Ceci est simple en utilisant le paquet R AlgDesign.

Pour la plupart des types de recherche, une telle conception serait suffisante. Cependant, dans l'étude de marché, on veut contrôler les effets d'ordre dans chaque bloc. C'est là que j'apprécierais de l'aide.

Créer des données de test

# The following code is not essential in understanding the problem, 
# but I provide it in case you are curious about the origin of the data itself. 
#library(AlgDesign) 
#set.seed(12345) 
#choices <- 4 
#nAttributes <- 7 
#blocksize <- 7 
#bsize <- rep(choices, blocksize) 
#PBIB <- optBlock(~., withinData=factor(1:nAttributes), blocksizes=bsize) 
#df <- data.frame(t(array(PBIB$rows, dim=c(choices, blocksize)))) 
#colnames(df) <- paste("Item", 1:choices, sep="") 
#rownames(df) <- paste("Set", 1:nAttributes, sep="") 

df <- structure(list(
    Item1 = c(1, 2, 1, 3, 1, 1, 2), 
    Item2 = c(4, 4, 2, 5, 3, 2, 3), 
    Item3 = c(5, 6, 5, 6, 4, 3, 4), 
    Item4 = c(7, 7, 6, 7, 6, 7, 5)), 
    .Names = c("Item1", "Item2", "Item3", "Item4"), 
    row.names = c("Set1", "Set2", "Set3", "Set4", "Set5", "Set6", "Set7"), 
    class = "data.frame") 

** Définir deux fonctions d'aide

balanceMatrix calcule l'équilibre de la matrice:

balanceMatrix <- function(x){ 
    t(sapply(unique(unlist(x)), function(i)colSums(x==i))) 
} 

balanceScore calcule une mesure de 'ajustement' - les scores les plus bas sont meilleurs, avec zéro parfait:

balanceScore <- function(x){ 
    sum((1-x)^2) 
} 

définir une fonction qui rééchantillonne les rangées au hasard

findBalance <- function(x, nrepeat=100){ 
    df <- x 
    minw <- Inf 
    for (n in 1:nrepeat){ 
     for (i in 1:nrow(x)){df[i,] <- sample(df[i, ])} 
     w <- balanceMatrix(df) 
     sumw <- balanceScore(w) 
     if(sumw < minw){ 
      dfbest <- df 
      minw <- sumw 
     } 
    } 
    dfbest 
} 

Main code

La trame de données df est une conception équilibrée de 7 ensembles. Chaque ensemble affichera 4 éléments au répondant. Les valeurs numériques dans df se réfère à 7 attributs différents. Par exemple, dans Set1, on demandera au répondant de choisir son option préférée parmi les attributs 1, 3, 4 et 7.

L'ordre des éléments dans chaque ensemble n'est conceptuellement pas important. Ainsi un ordre de (1,4,5,7) est identique à (7,5,4,1). Cependant, pour obtenir une conception entièrement équilibrée, chaque attribut apparaîtra un nombre égal de fois dans chaque colonne. Cette conception est là déséquilibrés, puisque l'attribut 1 apparaît 4 fois dans la colonne 1:

df 

    Item1 Item2 Item3 Item4 
Set1  1  4  5  7 
Set2  2  4  6  7 
Set3  1  2  5  6 
Set4  3  5  6  7 
Set5  1  3  4  6 
Set6  1  2  3  7 
Set7  2  3  4  5 

Pour essayer de trouver une conception plus équilibrée, je l'ai écrit la fonction findBalance. Cela conduit à une recherche aléatoire de meilleures solutions, en échantillonnant de manière aléatoire sur les lignes df. Avec 100 répétitions, il trouve la solution mieux suivante:

set.seed(12345) 
dfbest <- findBalance(df, nrepeat=100) 
dfbest 

    Item1 Item2 Item3 Item4 
Set1  7  5  1  4 
Set2  6  7  4  2 
Set3  2  1  5  6 
Set4  5  6  7  3 
Set5  3  1  6  4 
Set6  7  2  3  1 
Set7  4  3  2  5 

Cela semble plus équilibré, et la matrice de l'équilibre calculé contient beaucoup d'autres. La matrice de solde compte le nombre de fois que chaque attribut apparaît dans chaque colonne.Par exemple, le tableau suivant indique (dans la cellule supérieure gauche) cet attribut 1 apparaît deux fois pas du tout dans la colonne 1, et deux fois dans la colonne 2:

balanceMatrix(dfbest) 

    Item1 Item2 Item3 Item4 
[1,]  0  2  1  1 
[2,]  1  1  1  1 
[3,]  1  1  1  1 
[4,]  1  0  1  2 
[5,]  1  1  1  1 
[6,]  1  1  1  1 
[7,]  2  1  1  0 

Le solde score de pour cette solution est 6 , ce qui indique qu'il ya au moins six cellules inégales à 1:

balanceScore(balanceMatrix(dfbest)) 
[1] 6 

Ma question

Merci de suivre cet exemple détaillé. Ma question est comment puis-je réécrire cette fonction de recherche pour être plus systématique? Je voudrais dire R à:

• Réduire au minimum balanceScore(df)
• En changeant ligne ordre de df
• Sous réserve: déjà complètement contraint

Answer 1

OK, je en quelque sorte mal compris votre question. Alors, au revoir Fedorov, bonjour appliqué Fedorov.

L'algorithme suivant est basé sur la seconde itération de l'algorithme Fedorov:

1. calculer toutes de permutations possibles pour chaque ensemble, et de les stocker dans la liste de C0
2. tirer une première solution possible à partir du C0 espace (une permutation pour chaque ensemble). Cela peut être l'original, mais comme j'ai besoin des indices, je préfère commencer au hasard.
3. calculer le score pour chaque nouvelle solution, où le premier ensemble est remplacé par toutes les permutations.
4. remplacer le premier jeu avec la permutation donnant le score le plus bas
5. répétition 3-4 pour tous les autres ensemble
6. répétition 3-5 jusqu'à ce score atteint 0 ou pour n itérations.

Facultativement, vous pouvez redémarrer la procédure après 10 itérations et démarrer à partir d'un autre point de départ. En vous cas de test, il est apparu que quelques points de départ ont convergé très lentement à 0. La fonction ci-dessous trouvé des modèles expérimentaux équilibrés avec un score de 0 à en moyenne 1,5 secondes sur mon ordinateur:

> X <- findOptimalDesign(df) 
> balanceScore(balanceMatrix(X)) 
[1] 0 
> mean(replicate(20, system.time(X <- findOptimalDesign(df))[3])) 
[1] 1.733 

c'est donc le fonction maintenant (compte tenu de votre fonctions balanceMatrix et balanceScore d'origine):

findOptimalDesign <- function(x,iter=4,restart=T){ 
    stopifnot(require(combinat)) 
    # transform rows to list 
    sets <- unlist(apply(x,1,list),recursive=F) 
    nsets <- NROW(x) 
    # C0 contains all possible design points 
    C0 <- lapply(sets,permn) 
    n <- gamma(NCOL(x)+1) 

    # starting point 
    id <- sample(1:n,nsets) 
    Sol <- sapply(1:nsets,function(i)C0[[i]][id[i]]) 

    IT <- iter 
    # other iterations 
    while(IT > 0){ 
     for(i in 1:nsets){ 
      nn <- 1:n 
      scores <- sapply(nn,function(p){ 
      tmp <- Sol 
      tmp[[i]] <- C0[[i]][[p]] 
      w <- balanceMatrix(do.call(rbind,tmp)) 
      balanceScore(w) 
      }) 
      idnew <- nn[which.min(scores)] 
      Sol[[i]] <- C0[[i]][[idnew]] 

     } 
     #Check if score is 0 
     out <- as.data.frame(do.call(rbind,Sol)) 
     score <- balanceScore(balanceMatrix(out)) 
     if (score==0) {break} 
     IT <- IT - 1 

     # If asked, restart 
     if(IT==0 & restart){ 
      id <- sample(1:n,nsets) 
      Sol <- sapply(1:nsets,function(i)C0[[i]][id[i]]) 
      IT <- iter 
     } 
    } 
    out 
} 

HTH

EDIT: petit bug fixe (il redémarre immédiatement après chaque tour comme je l'ai oublié de condition sur IT). En faisant ça, ça marche encore un peu plus vite.