[R] library(SenoMineR)- Triangle Test Query

[Meyners, Michael]

Unless I missed it, neither the OP nor the list was CC'd on this, so for anyone interested, I forward this solution (untested from my side) from the package maintainer. Not sure whether the file comes through, so I include the updated code in the message's body below. 
Cheers, Michael


********** updated code for triangle.test ************
triangle.test <- function (design,answer,preference=NULL){

answer = gsub("(\\w)", "\\U\\1", as.character(answer), perl=TRUE)
labprod = levels(as.factor(c(as.character(design[,1]),as.character(design[,2]),as.character(design[,3]))))
nbprod = length(labprod)
nb.answer = nb.good = pref = matrix(0,nbprod,nbprod)
for (i in 1:nrow(design)){
  for (j in 1:nbprod){
     if (labprod[j] == design[i,1]) i1 = j
     if (labprod[j] == design[i,2]) i2 = j
     if (labprod[j] == design[i,3]) i3 = j
  }
  if (i1!=i2) nb.answer [i1,i2] = nb.answer[i1,i2]+1
  if (i1==i2) nb.answer [i1,i3] = nb.answer[i1,i3]+1
  if ((i1==i2)&(answer[i]=="Z")){
    nb.good[i1,i3]=nb.good[i1,i3]+1
    if (length(preference)>0){
      if (preference[i]!="Z") pref[i3,i1] = pref[i3,i1] +1
      if (preference[i]=="Z") pref[i1,i3] = pref[i1,i3] +1
    }
  }
  if ((i1==i3)&(answer[i]=="Y")){
    nb.good[i1,i2]=nb.good[i1,i2]+1
    if (length(preference)>0){
      if (preference[i]!="Y") pref[i2,i1] = pref[i2,i1] +1
      if (preference[i]=="Y") pref[i1,i2] = pref[i1,i2] +1
    }
  }
  if ((i2==i3)&(answer[i]=="X")){
    nb.good[i1,i2]=nb.good[i1,i2]+1
    if (length(preference)>0){
      if (preference[i]!="X") pref[i1,i2] = pref[i1,i2] +1
      if (preference[i]=="X") pref[i2,i1] = pref[i2,i1] +1
    }
  }
}
nb.good = nb.good + t(nb.good)
nb.answer = nb.answer + t(nb.answer)

diag(nb.answer)=1
prob = pbinom(nb.good-1,nb.answer,1/3,lower.tail=FALSE)
maxML = recognize = minimum = matrix(NA,nbprod,nbprod)
for (i in 1: (nbprod-1)){
  for (j in (i+1):nbprod){
    aux = matrix(0,nb.good[i,j]+1,1)
    for (k in 0:nb.good[i,j]) aux[k] = dbinom(nb.good[i,j]-k,nb.answer[i,j]-k,1/3)
    maxML[i,j] = maxML[j,i] = max(aux)
    recognize[i,j] = recognize[j,i] = rev(order(aux))[1]-1
    mini = 0
    for (k in round(nb.answer[i,j]/3):nb.answer[i,j]) if ((mini==0)&(dbinom(k,nb.answer[i,j],1/3)<0.05)) mini=k
    minimum[i,j]=minimum[j,i]=mini
  }
}

confusion = (nb.answer-recognize) / nb.answer
diag(nb.answer)=diag(recognize)=0
diag(maxML)=diag(confusion)=1

rownames(nb.answer) = colnames(nb.answer) = rownames(nb.good) = colnames(nb.good) = labprod
rownames(prob) = colnames(prob)= rownames(confusion) = colnames(confusion)= labprod
rownames(maxML) = colnames(maxML) = rownames(minimum) = colnames(minimum) = rownames(recognize) = colnames(recognize) = labprod
if (length(preference)>0) rownames(pref) = colnames(pref) = labprod

res = list()
res$nb.comp = nb.answer
res$nb.ident = nb.good
res$p.value = prob
res$nb.recognition = recognize
res$maxML = maxML
res$confusion = confusion
res$minimum = minimum
if (length(preference)>0) res$pref = pref
##res$complete = result
return(res)
}

********** end updated code for triangle.test ************




-----Original Message-----
From: Francois Husson 
Sent: Wednesday, June 08, 2011 14:43
To: Meyners, Michael
Subject: Re: [R] library(SenoMineR)- Triangle Test Query

  Dear Vijayan, dear Michael,

Indeed there was an error in the function triangle.test. I attach the new function.
Thanks Michael for your answer.
Best
Francois


Le 08/06/2011 12:32, Meyners, Michael a écrit :
> Vijayan,
>
> I cannot find an error in your code, but I had a look at the code of triangle.test -- unless I'm missing something, it contains a bug. If you study the way in which the matrix "pref" is updated, you find that the vector preference is compared to 1, 2 and 3 instead of "X", "Y" and "Z" as it should be. That way, some of the non-diagonal entries of pref will always be zero, irrespective of the data, which does not make sense. I think it should work if you modify the code accordingly. Alternatively, a quick patch (untested!) might be to code preferences as 1, 2 and 3 instead of the letters (but I'm not sure whether this has any other implications).
> I CC the author of the function and maintainer of the package; he should correct me if needed or could otherwise update the code for the next release (I worked on SensoMineR 1.11).
>
> Hope this helps,
> Michael
>
>> -----Original Message-----
>> From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-
>> project.org] On Behalf Of Vijayan Padmanabhan
>> Sent: Saturday, June 04, 2011 9:21
>> To: r-help at r-project.org
>> Subject: [R] library(SenoMineR)- Triangle Test Query
>>
>> Dear R Group
>> I was trying to use the triangle.test function in SensoMineR and
>> strangely i
>> encounter a error in the output of preference matrix from the analysis.
>> To illustrate, pl see the following dataframe of a design with the
>> response
>> and preference collected as shown below:
>>
>> design<-structure(list(`Product X` = c(3, 1, 4, 2, 4, 2, 1, 3, 4, 2,
>> 4, 2, 1, 3, 4, 2, 4, 2, 3, 1), `Product Y` = c(1, 1, 4, 4, 4,
>> 3, 1, 1, 4, 4, 4, 3, 1, 1, 4, 4, 4, 3, 1, 1), `Product Z` = c(3,
>> 2, 1, 2, 3, 3, 2, 3, 1, 2, 3, 3, 2, 3, 1, 2, 3, 3, 3, 2), Response =
>> structure(c(1L,
>> 2L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L,
>> 1L, 1L, 2L), .Label = c("X", "Z"), class = "factor"), Preference =
>> structure(c(1L,
>> 3L, 1L, 1L, 1L, 2L, 3L, 1L, 1L, 1L, 1L, 2L, 3L, 1L, 1L, 1L, 1L,
>> 2L, 1L, 2L), .Label = c("X", "Y", "Z"), class = "factor")), .Names =
>> c("Product X",
>> "Product Y", "Product Z", "Response", "Preference"), class =
>> "data.frame",
>> row.names = c("Panelist1.Test1",
>> "Panelist1.Test2", "Panelist2.Test1", "Panelist2.Test2",
>> "Panelist3.Test1",
>> "Panelist3.Test2", "Panelist4.Test1", "Panelist4.Test2",
>> "Panelist5.Test1",
>> "Panelist5.Test2", "Panelist6.Test1", "Panelist6.Test2",
>> "Panelist7.Test1",
>> "Panelist7.Test2", "Panelist8.Test1", "Panelist8.Test2",
>> "Panelist9.Test1",
>> "Panelist9.Test2", "Panelist10.Test1", "Panelist10.Test2"))
>>
>> If you were to investigate the above dataframe, you would find that for
>> the
>> comparision of Product 2 Vs Product 3, the preference indicates product
>> 3 is
>> preferred over product 2 all the time.
>>
>> ## Read output from the following script to see what i mean above:
>> subset(design,`Product X`==2&`Product Y`==3&`Product Z`==3)
>>
>> ##Output of above would be:
>> .                 Product X Product Y Product Z Response Preference
>> Panelist3.Test2         2         3         3        X          Y
>> Panelist6.Test2         2         3         3        X          Y
>> Panelist9.Test2         2         3         3        X          Y
>>
>> However when I analyse the design with the answers and preferences
>> using the
>> following script, I get the $pref output which shows that product 2 is
>> preferred over 3 all the time. Can somebody explain what is wrong in my
>> script?
>>
>> answer<-as.vector(design$Response)
>> preference<-as.vector(design$Preference)
>> triangle.test (design[,1:3], answer,preference)
>>
>> ##$pref output from the triangle.test function shows as follows:
>>
>> $pref
>>    1 2 3 4
>> 1 0 0 0 0
>> 2 4 0 3 0
>> 3 0 0 0 0
>> 4 0 0 0 0
>>
>>
>> Any help in helping me identify what is going wrong here would be
>> highly
>> appreciated.
>> Regards
>> Vijayan Padmanabhan
>>
>> 	[[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-
>> guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>