R/consistency_cata_panel.R

In ClustBlock: Clustering of Datasets

##=============================================================================


##' @title Test the consistency of the panel in a CATA experiment
##'
##' @description
##' Permutation test on the agreement between subjects in a CATA experiment
##'
##' @usage
##' consistency_cata_panel(Data,nblo, nperm=100, alpha=0.05)
##'
##'
##' @param Data data frame or matrix. Correspond to all the blocks of variables merged horizontally
##'
##' @param nblo  numerical. Number of blocks (subjects).
##'
##' @param nperm numerical. How many permutations are required? Default: 100
##'
##' @param alpha numerical between 0 and 1. What is the threshold? Default: 0.05
##'
##'
##'
##'@return a list with:
##'         \itemize{
##'          \item answer: the answer of the test
##'          \item pval: pvalue of the test
##'          \item dis: distance between the homogeneity and the median of the permutations
##'          }
##'
##'
##'
##' @keywords CATA RATA
##'
##' @references
##' Llobell, F., Giacalone, D., Labenne, A.,  Qannari, E.M. (2019).	Assessment of the agreement and cluster analysis of the respondents in a CATA experiment.	Food Quality and Preference, 77, 184-190.\cr
##' Bonnet, L., Ferney, T., Riedel, T., Qannari, E.M., Llobell, F. (September 14, 2022) .Using CATA for sensory profiling: assessment of the panel performance. Eurosense, Turku, Finland.
##'
##' @examples
##'\donttest{
##'  data(straw)
##' #with all subjects
##' consistency_cata_panel(Data=straw, nblo=114)
##'}
##'
##' @seealso   \code{\link{consistency_cata}}, \code{\link{change_cata_format}}, \code{\link{change_cata_format2}}
##'
##' @export


##=============================================================================
consistency_cata_panel=function(Data,nblo, nperm=100, alpha=0.05)
{
  nprod=nrow(Data)
  nattr=ncol(Data)/nblo
  n=nprod
  Blocks=rep(nattr,nblo)
  J=rep(1:nblo , times =  Blocks )# indicates which block each variable belongs to
  X=Data

  #parapet for numerical Data
  for (i in 1: ncol(Data))
  {
    if (is.numeric(Data[,i])==FALSE)
    {
      stop(paste("The data must be numeric (column",i,")"))
    }
  }

  # #parapet for binary Data
  # if ((sum(Data==0)+sum(Data==1))!=(dim(Data)[1]*dim(Data)[2]))
  # {
  #   stop("only binary Data is accepted (0 or 1)")
  # }

  #no NA
  if(sum(is.na(Data))>0)
  {
    print("NA detected:")
    tabna=which(is.na(Data), arr.ind = TRUE)
    print(tabna)
    stop(paste("NA are not accepted"))
  }


  Xj=array(0,dim=c(n,nattr,nblo))  # array with all subjects matrices
  muk=NULL
  for(j in 1:nblo)
  {
    Aj=as.matrix(X[,J==j])
    normXj=sqrt(sum(diag(tcrossprod(Aj, Aj))))
    muk[j]=normXj
    if(normXj==0)
    {
      stop(paste("error: the subject", j, "has only 0 or only 1"))  #parapet for constant configurations
    }
    Xj[,,j]=Aj/normXj #standardization
  }


  # S matrix:
  S=matrix(0,nblo,nblo)
  diag(S)=rep(1,nblo)
  for (i in 1:(nblo-1)) {
    for (j in (i+1):nblo) {
      S[i,j]=sum(diag(tcrossprod(Xj[,,i],Xj[,,j])))
      S[j,i]=S[i,j]
    } }


  #first eigenvector and eigenvalue of S matrix
  ressvd=svd(S)
  lambda=ressvd$d[1]
  hom=lambda/nblo

  #permutations
  hom_perm=NULL
  lambda_perm=NULL
  for (k in 1:nperm)
  {

    for(j in 1:nblo)
    {
      Aj_perm=Aj=as.matrix(Xj[,,j])
      tirage=sample(1:nrow(X))
      for (i in 1:nrow(X))
      {
        Aj_perm[i,]=Aj[tirage[i],]
      }

      Xj[,,j]=Aj_perm
    }


    # S matrix:
    S=matrix(0,nblo,nblo)
    diag(S)=rep(1,nblo)
    for (i in 1:(nblo-1)) {
      for (j in (i+1):nblo) {
        S[i,j]=sum(diag(tcrossprod(Xj[,,i],Xj[,,j])))
        S[j,i]=S[i,j]
      } }



    #first eigenvector and eigenvalue of S matrix
    ressvd=svd(S)
    lambda_perm[k]=ressvd$d[1]
    hom_perm[k]=lambda_perm[k]/nblo
  }
  pval=sum(hom_perm>=hom)/nperm
  if (pval<alpha)
  {
    consist="the panel is consistent"
  }else{
    consist="the panel is not consistent"
  }

  dis=hom-median(hom_perm)

  return(list(answer= consist, pval=pval, dis=dis))
}