Problematic time interpolation for user defined PK parameters

[Yuri05]

For the user defined PK parameters with a given time range, the indices of the start/end value points are estimated using the ArrayHelper.ClosestIndexOf() function :

OSPSuite.Core/src/OSPSuite.Core/Domain/Services/PKValuesCalculator.cs

Lines 159 to 160 in fce3308

	var startIndex = ArrayHelper.ClosestIndexOf(time, startTime);
	var endIndex = ArrayHelper.ClosestIndexOf(time, endTime);

(also used for the estimation of time indices for the dosing intervals)

OSPSuite.Core/src/OSPSuite.Core/Domain/Services/PKValuesCalculator.cs

Lines 197 to 198 in fce3308

	var dosingStartIndex = ArrayHelper.ClosestIndexOf(time, dosingInterval.StartValue);
	var dosingEndIndex = ArrayHelper.ClosestIndexOf(time, dosingInterval.EndValue);

The function ArrayHelper.ClosestIndexOf()is defined as following:

OSPSuite.Core/src/OSPSuite.Core/Domain/Services/PKValuesCalculator.cs

Lines 267 to 278 in fce3308

	private static int closestIndexOf(List<float> list, float value)
	{
	int index = list.BinarySearch(value);
	//already exist return
	if (index >= 0)
	return index;
	//does not exist in the list. We find the index of the point immediately after this value
	var closestIndex = ~index;
	return closestIndex < list.Count ? closestIndex : -1;
	}
	}

The problem is: if the time point is not on the simulation output raster, the function returns not the closest time index but the next time index after the time of interest.

This might result in significant deviations compared to the expected result.

E.g. in the attached simulation, the output schema results in the time vector (in days):

...
12,99848
13,00889
...

46_Ctrough.zip

The applications are given once daily (particularly at t=13 days).

Simulated values for Peripheral Venous Blood look like this:

Time [day(s)]	PeripheralVenousBlood-C1-Plasma (Peripheral Venous Blood) [µmol/l]
…	…
12,99848	1,624205
13,00889	3,23945
…	…

Now when defining a user defined PK parameter for CTrough between 12 and 13 days, I would expect to get the value 1,62. But I get 3,24 instead.

library(ospsuite)

sim <- loadSimulation("S1.pkml")
simRes <- runSimulation(sim)

ctrough_md <- addUserDefinedPKParameter(name = "ctrough_md", standardPKParameter = StandardPKParameter$C_trough)
ctrough_md$startTime <- 12*24*60 #17280 min
ctrough_md$endTime <- 13*24*60   #18720 min

pkAnalysis <- calculatePKAnalyses(results = simRes)
outputPath <- simRes$allQuantityPaths[[1]]
print(outputPath)
allPkParams <- pkAnalysis$allPKParametersFor(outputPath)

trough <- pkAnalysis$pKParameterFor(quantityPath = outputPath, pkParameter = "ctrough_md")
trough$values

[PavelBal]

I think that the behavior of ArrayHelper.ClosestIndexOf() is correct, as the distance from 13 to 13,00889 is smaller than from 13 to 12,99848. The PKValuesCalculator should do some additional magic, e.g., (in R pseudo code):

if  (time[[endIndex]] > endTime) {
  # Return the value that is within the defined interval
    return(endIndex <- endIndex - 1)
}
# If the closest time is lower than the start time, than the results resolution 
# is not sufficient to calculate the value. Error?
if (time[[endIndex]] < startTime) {
  error()
}
return(endIndex)

[Yuri05]

the distance from 13 to 13,00889 is smaller than from 13 to 12,99848.

ahm, not really 😉

Apart from that: after discussion with @AnnikaRPS the suggestion is to create a Warning if start time and/or end time of a user defined PK parameter does not belong to the output time raster of the simulation. In 99% of all cases it does not really matter, but in few cases like the one in the issue description the deviation is really significant (which could happen even if the correct index is chosen or the 2 adjacent) simulated values are interpolated)

[AnnikaRPS]

I agree with @Yuri05 , that a warning is needed.
Also, @PavelBal, I don't think the function is doing the correct thing since the user explicitly gave the end time of 13 days, thus excluding all timepoints >13days. This should have a higher priority than e.g. 13.1 being closer to 13 than 12.7. Time values >13 should be excluded due to the user defined filter.

[PavelBal]

ahm, not really

Ah agree, sorry I misread the numbers ;)

But still

I don't think the function is doing the correct thing since the user explicitly gave the end time of 13 days, thus excluding all timepoints >13days.

I am talking about the function closestIndexOf() which does not include any boundaries. But obviously its implementation is wrong.

btw this is how we implemented this in R:

getIndexClosestToValue <- function(value, array, thresholdAbs = NULL, thresholdRel = NULL) {
  # If no absolute threshold is set, calculate if from relative threshold
  if (is.null(thresholdAbs)) {
    # If no relative threshold is set also, no threshold is applied
    if (!is.null(thresholdRel)) {
      thresholdAbs <- abs(value * thresholdRel)
    } else {
      thresholdAbs <- Inf
    }
  }

  validateIsNumeric(c(value, array))

  # Calculate distances
  distances <- abs(array - value)
  idx <- which(distances == min(distances) & distances <= thresholdAbs)

  if (length(idx) == 0) {
    warning(messages$warningValueWithinThresholdNotExisting(value, thresholdAbs))
    return(NULL)
  }

  return(idx)
}