merge to master

R/Dockerfile

@@ -9,13 +9,15 @@ RUN DEBIAN_FRONTEND=noninteractive apt install -y -q libgmp-dev
 
 RUN R -e 'install.packages("hitandrun", repos="http://cran.rstudio.com/"); if (!require("hitandrun")) quit(save="no", status=8)'
 RUN R -e 'install.packages("smaa", repos="http://cran.rstudio.com/"); if (!require("smaa")) quit(save="no", status=8)'
+RUN R -e 'install.packages("abind", repos="http://cran.rstudio.com/"); if (!require("abind")) quit(save="no", status=8)'
 
 ADD *.R /tmp/
 ADD util/*.R /tmp/
+RUN rm  /tmp/apiEntryPoint.R
 
 RUN cat /tmp/*.R > /var/lib/patavi/smaa_service.R
 
 USER patavi
 WORKDIR /var/lib/patavi
 
-ENTRYPOINT ["patavi-worker", "--method", "smaa_v2", "-n", "1", "--file", "/var/lib/patavi/smaa_service.R", "--rserve", "--packages", "MASS,rcdd,hitandrun,smaa"]
+ENTRYPOINT ["patavi-worker", "--method", "smaa_v2", "-n", "1", "--file", "/var/lib/patavi/smaa_service.R", "--rserve", "--packages", "MASS,rcdd,hitandrun,smaa,abind"]

R/apiEntryPoint.R

@@ -4,6 +4,7 @@ library(hitandrun)
 
 smaa_v2 <- function(params) {
   allowed <- c(
+    'choiceBasedMatching',
     'deterministic',
     'indifferenceCurve',
     'matchingElicitationCurve',

R/choiceBasedMatching.R

@@ -0,0 +1,250 @@
+run_choiceBasedMatching <- function(params) {
+  numberOfCriteria <- length(params[['criteria']])
+  criterionIds <- c()
+  for (criterion in params[['criteria']]) {
+    criterionIds <- append(criterionIds, criterion[['id']])
+  }
+
+  pvf <- lapply(params[['criteria']], createPvf)
+  names(pvf) <- criterionIds
+
+  constraintsFromHistory <- generateConstraintsFromHistory(params, criterionIds, pvf)
+  edgeLengths <- calculateEdgeLengths(constraintsFromHistory)
+
+  fromR <- params
+
+  if (isDoneEliciting(params)) {
+    fromR[['preferences']] <- retrieveUpperBoundConstraints(constraintsFromHistory, numberOfCriteria, criterionIds)
+  } else {
+    nextQuestion <- generateNextQuestion(params, constraintsFromHistory, edgeLengths, criterionIds, pvf)
+    fromR[['answersAndQuestions']] <- append(fromR[['answersAndQuestions']], list(list(question = nextQuestion)))
+  }
+  return(fromR)
+}
+
+isDoneEliciting <- function(params) {
+  numberOfCriteria <- length(params[['criteria']])
+  numberOfAnswers <- length(params[['answersAndQuestions']])
+  if (numberOfCriteria <= 4) {
+    return(numberOfAnswers == (numberOfCriteria - 1) * 4)
+  } else {
+    return(numberOfAnswers == (numberOfCriteria - 2) * 4)
+  }
+}
+
+retrieveUpperBoundConstraints <- function(constraint, numberOfCriteria, criterionIds) {
+  numberOfUpperBoundConstraints <- nrow(constraint[['constr']]) - (numberOfCriteria + 1)
+  preferences <- vector("list", numberOfUpperBoundConstraints)
+  for (index in 1:numberOfUpperBoundConstraints) {
+    criteria <- getCriteriaWithConstraints(criterionIds, constraint, index, numberOfCriteria)
+    bound <- getBound(constraint, index, numberOfCriteria)
+    preferences[[index]] <- list(type = "upper ratio", elicitationMethod = "choice", criteria = criteria, bound = bound)
+  }
+  return(preferences)
+}
+
+getCriteriaWithConstraints <- function(criterionIds, constraint, index, numberOfCriteria) {
+  criterion1Id <- criterionIds[getConstraint(constraint, index, numberOfCriteria) > 0]
+  criterion2Id <- criterionIds[getConstraint(constraint, index, numberOfCriteria) < 0]
+  return(c(criterion1Id, criterion2Id))
+}
+
+getBound <- function(constraint, index, numberOfCriteria) {
+  bound <- -1 * constraint[['constr']][index + (numberOfCriteria + 1), getConstraint(constraint, index, numberOfCriteria) < 0] /
+    constraint[['constr']][index + (numberOfCriteria + 1), getConstraint(constraint, index, numberOfCriteria) > 0]
+  return(bound)
+}
+
+getConstraint <- function(constraint, index, numberOfCriteria) {
+  return(constraint[['constr']][index + (numberOfCriteria + 1),])
+}
+
+generateNextQuestion <- function(params, constraints, edgeLengths, criterionIds, pvf) {
+  criterionIdsEdgeIndices <- calculateEdgeIndices(edgeLengths)
+  criterionIdsEdge <- criterionIds[criterionIdsEdgeIndices]
+  cutPoint <- calculateCutPoint(constraints, criterionIdsEdgeIndices)
+
+  alternativeA <- rep(0, 2)
+  names(alternativeA) <- c("criterion1Value", "criterion2Value")
+  alternativeB <- alternativeA
+
+  criterion1Range <- params[['criteria']][[which(criterionIds == criterionIdsEdge[1])]][['pvf']][['range']]
+  criterion2Range <- params[['criteria']][[which(criterionIds == criterionIdsEdge[2])]][['pvf']][['range']]
+
+  alternativeA['criterion1Value'] <- selectQuestionValueFromRange(params, criterionIds, criterionIdsEdge[1], criterion1Range, 1, 2)
+  alternativeB['criterion1Value'] <- selectQuestionValueFromRange(params, criterionIds, criterionIdsEdge[1], criterion1Range, 2, 1)
+
+  alternativeA['criterion2Value'] <- selectQuestionValueFromRange(params, criterionIds, criterionIdsEdge[2], criterion2Range, 2, 1)
+  alternativeB['criterion2Value'] <- selectQuestionValueFromRange(params, criterionIds, criterionIdsEdge[2], criterion2Range, 1, 2)
+
+  if (cutPoint >= 1) {
+    alternativeB['criterion1Value'] <- calculateQuestionValue(pvf, criterionIdsEdge[1], 1 / cutPoint, criterion1Range)
+  } else {
+    alternativeA['criterion2Value'] <- calculateQuestionValue(pvf, criterionIdsEdge[2], cutPoint, criterion2Range)
+  }
+
+  isFirstAlternativeA <- runif(1) <= 0.5
+  if (isFirstAlternativeA) {
+    question <- list(A = alternativeA, B = alternativeB, criterionIds = criterionIdsEdge)
+  } else {
+    question <- list(A = alternativeB, B = alternativeA, criterionIds = criterionIdsEdge)
+  }
+  return(question)
+}
+
+calculateEdgeIndices <- function(edgeLengths) {
+  candidateEdges <- which(edgeLengths == max(edgeLengths), arr.ind = T)
+  edge <- sample.int(nrow(candidateEdges), 1)
+  criterionIdsEdgeIndices <- c(min(candidateEdges[edge,]), max(candidateEdges[edge,]))
+  return(criterionIdsEdgeIndices)
+}
+
+calculateCutPoint <- function(constraints, criterionIdsEdgeIndices) {
+  samples <- hitandrun(constraints, 1e2)
+  cutPoint <- median(samples[, criterionIdsEdgeIndices[1]] / samples[, criterionIdsEdgeIndices[2]])
+  return(cutPoint)
+}
+
+selectQuestionValueFromRange <- function(params, criterionIds, criterionId, criterionRange, index1, index2) {
+  if (params[['criteria']][[which(criterionIds == criterionId)]][['pvf']][['direction']] == "increasing") {
+    value <- criterionRange[index1]
+  } else {
+    value <- criterionRange[index2]
+  }
+  return(value)
+}
+
+calculateQuestionValue <- function(pvf, criterionId, cutPointValue, criterionRange) {
+  value <- uniroot(f = function(x) { pvf[[criterionId]](x) - cutPointValue }, interval = criterionRange)[['root']]
+  stepSize <- 10 ^ floor(log10(diff(criterionRange)) - 1)
+  roundedValue <- round(value / stepSize, 0) * stepSize
+  return(roundedValue)
+}
+
+generateConstraintsFromHistory <- function(params, criterionIds, pvf) {
+  numberOfCriteria <- length(params[['criteria']])
+  constraints <- simplexConstraints(numberOfCriteria)
+  if (length(params[['answersAndQuestions']]) > 0) {
+    for (answerAndQuestion in params[['answersAndQuestions']]) {
+      currentQuestionConstraint <- generateCurrentQuestionConstraint(answerAndQuestion, numberOfCriteria, criterionIds, pvf)
+      constraints <- mergeConstraints(constraints, currentQuestionConstraint)
+    }
+  }
+  return(constraints)
+}
+
+generateCurrentQuestionConstraint <- function(answerAndQuestion, numberOfCriteria, criterionIds, pvf) {
+  currentQuestionConstraint <- rep(0, numberOfCriteria)
+
+  criterionId1 <- which(answerAndQuestion[['question']][['criterionIds']][1] == criterionIds)
+  currentQuestionConstraint[criterionId1] <- calculateConstraint(pvf, criterionId1, answerAndQuestion, "criterion1Value")
+
+  criterionId2 <- which(answerAndQuestion[['question']][['criterionIds']][2] == criterionIds)
+  currentQuestionConstraint[criterionId2] <- calculateConstraint(pvf, criterionId2, answerAndQuestion, "criterion2Value")
+
+  if (answerAndQuestion[['answer']] == "B") {
+    currentQuestionConstraint <- -1 * currentQuestionConstraint
+  }
+  currentQuestionConstraint <- list(constr = currentQuestionConstraint, dir = "<=", rhs = 0)
+  return(currentQuestionConstraint)
+}
+
+calculateConstraint <- function(pvf, criterionId, answerAndQuestion, criterionSelector) {
+  constraint <- pvf[[criterionId]](answerAndQuestion[['question']][['B']][criterionSelector]) - pvf[[criterionId]](answerAndQuestion[['question']][['A']][criterionSelector])
+  return(constraint)
+}
+
+euclideanDistance <- function(x, y) {
+  stopifnot(length(x) == length(y))
+  result <- sqrt(sum((x - y) ^ 2))
+  return(result)
+}
+
+caculateRatioBounds <- function(constr, index1, index2) {
+  # Guide on how to transform the constraint set to maximize or minimize a weight ratio wi/wj: http://lpsolve.sourceforge.net/5.1/ratio.htm 
+  A <- cbind(-constr[['rhs']], constr[['constr']])
+  b <- rep(0, length(constr[['rhs']]))
+  tranformedConstraint <- list(constr = A, rhs = b, dir = constr[['dir']])
+
+  c <- rep(0, ncol(constr[['constr']]))
+  c[index2] <- 1
+  cTransformed <- c(0, c)
+  y0Constr <- list(constr = cTransformed, rhs = 1, dir = "=")
+
+  tranformedConstraint <- mergeConstraints(tranformedConstraint, y0Constr)
+  hrepTransformed <- makeH(a1 = tranformedConstraint[['constr']][tranformedConstraint[['dir']] == "<=",], b1 = tranformedConstraint[['rhs']][tranformedConstraint[['dir']] == "<="],
+                            a2 = tranformedConstraint[['constr']][tranformedConstraint[['dir']] == "=",], b2 = tranformedConstraint[['rhs']][tranformedConstraint[['dir']] == "="])
+
+  objectiveFunction <- rep(0, length(cTransformed))
+  objectiveFunction[index1 + 1] <- 1
+
+  # Obtain upper bound for wi/wj
+  upper <- lpcdd(hrepTransformed, objectiveFunction, minimize = F)
+  if (upper[['solution.type']] == "Optimal") {
+    weightsUpper <- upper[['primal.solution']][2:length(upper[['primal.solution']])] / upper[['primal.solution']][1]
+    upperRatioBound <- weightsUpper[index1] / weightsUpper[index2]
+  } else {
+    upperRatioBound <- Inf
+  }
+
+  # Obtain lower bound for wi/wj
+  lower <- lpcdd(hrepTransformed, objectiveFunction, minimize = T)
+  weightsLower <- lower[['primal.solution']][2:length(lower[['primal.solution']])] / lower[['primal.solution']][1]
+  lowerRatioBound <- weightsLower[index1] / weightsLower[index2]
+
+  ratioBounds <- c(lowerRatioBound, upperRatioBound)
+  return(ratioBounds)
+}
+
+obtainAllRatioBounds <- function(constraints) {
+  numberOfAttributes <- ncol(constraints[['constr']])
+  ratioBounds <- array(dim = c(numberOfAttributes, numberOfAttributes, 2))
+  for (index1 in 1:numberOfAttributes) {
+    for (index2 in 1:numberOfAttributes) {
+      if (index1 == index2) {
+        ratioBounds[index1, index2,] <- c(1, 1)
+      } else {
+        ratioBounds[index1, index2,] <- caculateRatioBounds(constraints, index1, index2)
+      }
+    }
+  }
+  return(ratioBounds)
+}
+
+calculateEdgeRatioBoundIntersection <- function(numberOfcriteria, index1, index2, bound) {
+  # Hyperplane: wi / wj = bound
+  # Formulas taken from: https://en.wikipedia.org/wiki/Line%E2%80%93plane_intersection
+
+  pointHyperplane <- rep(0, numberOfcriteria) # origin is always on plane
+  vertex1 <- pointHyperplane
+  vertex1[index1] <- 1
+  vertex2 <- pointHyperplane
+  vertex2[index2] <- 1
+
+  if (bound == 0) {
+    intersectionPoint <- vertex2
+  } else {
+    if (bound == Inf) {
+      intersectionPoint <- vertex1
+    } else {
+      normalHyperplane <- rep(0, numberOfcriteria)
+      normalHyperplane[index1] <- 1
+      normalHyperplane[index2] <- -bound
+      intersectionPoint <- vertex1 + as.numeric(((pointHyperplane - vertex1) %*% normalHyperplane) / ((vertex1 - vertex2) %*% normalHyperplane)) * (vertex1 - vertex2)
+    }
+  }
+  return(intersectionPoint)
+}
+
+calculateEdgeLengths <- function(constraints) {
+  numberOfcriteria <- ncol(constraints[['constr']])
+  edgeLengths <- matrix(0, nrow = numberOfcriteria, ncol = numberOfcriteria)
+  allRatioBounds <- obtainAllRatioBounds(constraints)
+  for (index1 in 1:(numberOfcriteria - 1)) {
+    for (index2 in (index1 + 1):numberOfcriteria) {
+      edgeLengths[index1, index2] <- euclideanDistance(calculateEdgeRatioBoundIntersection(numberOfcriteria, index1, index2, allRatioBounds[index1, index2, 1]), calculateEdgeRatioBoundIntersection(numberOfcriteria, index1, index2, allRatioBounds[index1, index2, 2]))
+      edgeLengths[index2, index1] <- edgeLengths[index1, index2]
+    }
+  }
+  return(round(edgeLengths, 3))
+}

R/smaa.R

@@ -53,7 +53,7 @@ applyMeasurementUncertainty <- function(params, criteria, measurements) {
   for (criterion in criteria) {
     medianMeasurements[, criterion] <- pvf[[criterion]](medianMeasurements[, criterion])
   }
-  for (i in 1:hitAndRunSamples) {
+  for (i in 1:dim(measurements)[1]) {
     measurements[i,,] <- medianMeasurements
   }
   return(measurements)

R/util/constraint.R

@@ -9,6 +9,8 @@ genHARconstraint <- function(statement, criteria) {
     return(getRatioBoundConstraint(statement$bounds, numberOfCriteria, index1, index2))
   } else if (statement$type == "exact swing") {
     return(getRatioConstraint(numberOfCriteria, index1, index2, statement$ratio))
+  } else if (statement$type == "upper ratio") {
+    return(upperRatioConstraint(numberOfCriteria,  index1,  index2, statement$bound))
   }
 }
 

R/util/pvf.R

@@ -1,5 +1,4 @@
 # import sample form sampler.R
-# import createPvf from pvf.R
 
 createPvf <- function(criterion) {
   pvf <- criterion$pvf

app/ts/McdaApp/Workspace/CurrentScenarioContext/CurrentScenarioContext.tsx

@@ -115,7 +115,10 @@ export function CurrentScenarioContextProviderComponent({
 
   const getWeights = useCallback(
     (scenario: IMcdaScenario, pvfs: Record<string, TPvf>): void => {
-      if (scenario.state.prefs[0]?.elicitationMethod === 'imprecise') {
+      if (
+        scenario.state.prefs[0]?.elicitationMethod === 'imprecise' ||
+        scenario.state.prefs[0]?.elicitationMethod === 'choice'
+      ) {
         getWeightsFromPatavi(scenario, pvfs);
       } else {
         getWeightsThroughCalculation(scenario);

app/ts/McdaApp/Workspace/CurrentTab/Preferences/Preferences.tsx

@@ -15,14 +15,17 @@ import {EquivalentChangeContextProviderComponent} from './EquivalentChange/Equiv
 import AdvancedPartialValueFunction from './PartialValueFunctions/AdvancedPartialValueFunctions/AdvancedPartialValueFunction';
 import {AdvancedPartialValueFunctionContextProviderComponent} from './PartialValueFunctions/AdvancedPartialValueFunctions/AdvancedPartialValueFunctionContext/AdvancedPartialValueFunctionContext';
 import PreferencesView from './PreferencesView/PreferencesView';
+import {ErrorContext} from 'app/ts/Error/ErrorContext';
 
 export default function Preferences() {
   const {filteredCriteria, stepSizesByCriterion} = useContext(
     CurrentSubproblemContext
   );
   const {setActiveView, currentScenario, activeView, pvfs, updateScenario} =
     useContext(CurrentScenarioContext);
-  const {showPercentages} = useContext(SettingsContext);
+  const {showPercentages, showCbmPieChart} = useContext(SettingsContext);
+  const {setErrorMessage} = useContext(ErrorContext);
+
   const {
     workspace: {
       properties: {title}
@@ -64,6 +67,9 @@ export default function Preferences() {
       case 'ranking':
         document.title = 'Ranking';
         break;
+      case 'choice':
+        document.title = 'Choice-based matching';
+        break;
       case 'threshold':
         document.title = 'Threshold technique elicitation';
         break;
@@ -98,6 +104,8 @@ export default function Preferences() {
           manualLexicon={lexicon}
           manualHost={'@MCDA_HOST'}
           manualPath="/manual.html"
+          showCbmPieChart={showCbmPieChart}
+          setErrorMessage={setErrorMessage}
         />
       );
     } else {

app/ts/McdaApp/Workspace/CurrentTab/Preferences/PreferencesWeights/PreferencesWeightsButtons/PreferencesWeightsButtons.tsx

@@ -32,6 +32,9 @@ export default function PreferencesWeightsButtons() {
     setActiveView('imprecise');
   }
 
+  function handleChoiceBasedClick() {
+    setActiveView('choice');
+  }
   function handleThresholdClick() {
     setActiveView('threshold');
   }
@@ -88,6 +91,22 @@ export default function PreferencesWeightsButtons() {
           Imprecise Swing Weighting
         </Button>
       </Tooltip>
+      <Tooltip
+        title={
+          containsNonLinearPvf
+            ? ''
+            : 'Saving this preference will reset all criteria trade-off preferences'
+        }
+      >
+        <Button
+          id="choice-based-matching-button"
+          onClick={handleChoiceBasedClick}
+          color="primary"
+          variant="contained"
+        >
+          Choice-based Matching
+        </Button>
+      </Tooltip>
       <Tooltip
         title={
           containsNonLinearPvf

app/ts/McdaApp/Workspace/ScenariosContext/TElicitationMethod.ts

@@ -3,5 +3,6 @@ export type TElicitationMethod =
   | 'precise'
   | 'matching'
   | 'imprecise'
+  | 'choice'
   | 'threshold'
   | 'none';

app/ts/McdaApp/Workspace/ScenariosContext/TPreferencesView.ts

@@ -4,5 +4,6 @@ export type TPreferencesView =
   | 'precise'
   | 'imprecise'
   | 'matching'
-  | 'advancedPvf'
-  | 'threshold';
+  | 'choice'
+  | 'threshold'
+  | 'advancedPvf';