Source code for pkpd.protocols.protocol_pkpd_ode_bootstrap

# **************************************************************************
# *
# * Authors:     Carlos Oscar Sorzano (
# *
# * Kinestat Pharma
# *
# * This program is free software; you can redistribute it and/or modify
# * it under the terms of the GNU General Public License as published by
# * the Free Software Foundation; either version 2 of the License, or
# * (at your option) any later version.
# *
# * This program is distributed in the hope that it will be useful,
# * but WITHOUT ANY WARRANTY; without even the implied warranty of
# * GNU General Public License for more details.
# *
# * You should have received a copy of the GNU General Public License
# * along with this program; if not, write to the Free Software
# * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
# * 02111-1307  USA
# *
# *  All comments concerning this program package may be sent to the
# *  e-mail address ''
# *
# **************************************************************************

import numpy as np

import pyworkflow.protocol.params as params
from pkpd.objects import PKPDFitting, PKPDSampleFitBootstrap, PKPDLSOptimizer
from pyworkflow.protocol.constants import LEVEL_ADVANCED
from .protocol_pkpd_ode_base import ProtPKPDODEBase


[docs]class ProtPKPDODEBootstrap(ProtPKPDODEBase): """ Bootstrap of an ODE protocol""" _label = 'ODE bootstrap' PRM_PROTOCOL = 0 PRM_FITTING = 1 #--------------------------- DEFINE param functions -------------------------------------------- def _defineParams(self, form): form.addSection('Input') form.addParam('inputODE', params.PointerParam, label="Input ODE model", pointerClass='ProtPKPDMonoCompartment, ProtPKPDMonoCompartmentUrine, ProtPKPDTwoCompartments, '\ 'ProtPKPDTwoCompartmentsAutoinduction, ProtPKPDTwoCompartmentsClint, '\ 'ProtPKPDTwoCompartmentsClintCl, '\ 'ProtPKPDTwoCompartmentsClintMetabolite, ProtPKPDTwoCompartmentsUrine, '\ 'ProtPKPDODERefine', help='Select a run of an ODE model') form.addParam('paramsSource', params.EnumParam, label="Source of parameters", choices=['ODE protocol','ODE Fitting'], default=0, help="Choose a population of parameters, your own parameters or a previously fitted set of measurements") form.addParam('inputFitting', params.PointerParam, label="Input ODE fitting", condition="paramsSource==1", pointerClass='PKPDFitting', help='It must be a fitting coming from a compartmental PK fitting') form.addParam('Nbootstrap', params.IntParam, label="Bootstrap samples", default=200, expertLevel=LEVEL_ADVANCED, help='Number of bootstrap realizations for each sample') form.addParam('sampleLength', params.IntParam, label="Sample length", default=-1, expertLevel=LEVEL_ADVANCED, help='If the input experiment represents a population, the bootstrap sample is so overdetermined '\ 'that the bootstrap parameter estimate seldom moves from the original values. In this case, '\ 'you may use the sample length to generate bootstrap samples of the same length as the indiviudals '\ 'taking part of the population. If this value is set to -1, then the length of the bootstrap sample '\ 'will be the same as the one in the input experiment. If set to any other value, e.g. 8, then '\ 'each bootstrap sample will have this length') form.addParam('confidenceInterval', params.FloatParam, label="Confidence interval", default=95, expertLevel=LEVEL_ADVANCED, help='Confidence interval for the fitted parameters') form.addParam('deltaT', params.FloatParam, default=2, label='Step (min)', expertLevel=LEVEL_ADVANCED) #--------------------------- INSERT steps functions -------------------------------------------- def _insertAllSteps(self): self._insertFunctionStep('runFit',self.inputODE.get().getObjId(), self.Nbootstrap.get(), self.confidenceInterval.get()) self._insertFunctionStep('createOutputStep') #--------------------------- STEPS functions --------------------------------------------
[docs] def parseBounds(self, boundsString): if boundsString!="" and boundsString!=None: tokens=boundsString.split(';') if len(tokens)!=self.getNumberOfParameters(): raise Exception("The number of bound intervals does not match the number of parameters") self.boundsList=[] for token in tokens: values = token.strip().split(',') self.boundsList.append((float(values[0][1:]),float(values[1][:-1])))
[docs] def setBounds(self,sample): self.parseBounds(self.protODE.bounds.get()) self.setBoundsFromBoundsList()
[docs] def setBoundsFromBoundsList(self): Nbounds = len(self.boundsList) Nsource = self.drugSource.getNumberOfParameters() Nmodel = self.model.getNumberOfParameters() if Nbounds!=Nsource+Nmodel: raise "The number of parameters (%d) and bounds (%d) are different"%(Nsource+Nmodel,Nbounds) self.boundsSource = self.boundsList[0:Nsource] self.boundsPK = self.boundsList[Nsource:] self.model.bounds = self.boundsPK
[docs] def getBounds(self): return self.boundsList
[docs] def runFit(self, objId, Nbootstrap, confidenceInterval): self.protODE = self.inputODE.get() if self.paramsSource.get()==self.PRM_PROTOCOL: self.experiment = self.readExperiment(self.protODE.outputExperiment.fnPKPD) self.fitting = self.readFitting(self.protODE.outputFitting.fnFitting) else: self.fitting = self.readFitting(self.inputFitting.get().fnFitting) self.experiment = self.readExperiment(self.inputFitting.get().fnExperiment) # Get the X and Y variable names self.varNameX = self.fitting.predictor.varName if type(self.fitting.predicted)==list: self.varNameY = [v.varName for v in self.fitting.predicted] else: self.varNameY = self.fitting.predicted.varName self.protODE.experiment = self.experiment self.protODE.varNameX = self.varNameX self.protODE.varNameY = self.varNameY # Create output object self.fitting = PKPDFitting("PKPDSampleFitBootstrap") self.fitting.fnExperiment.set(self.experiment.fnPKPD.get()) self.fitting.predictor=self.experiment.variables[self.varNameX] if type(self.varNameY)==list: self.fitting.predicted=[self.experiment.variables[v] for v in self.varNameY] else: self.fitting.predicted=self.experiment.variables[self.varNameY] self.fitting.modelParameterUnits = None # Actual fitting if self.protODE.fitType.get()==0: fitType = "linear" elif self.protODE.fitType.get()==1: fitType = "log" elif self.protODE.fitType.get()==2: fitType = "relative" parameterNames = None for groupName, group in self.experiment.groups.items(): self.printSection("Fitting "+groupName) self.protODE.clearGroupParameters() self.clearGroupParameters() for sampleName in group.sampleList: print(" Sample "+sampleName) sample = self.experiment.samples[sampleName] self.protODE.createDrugSource() self.protODE.setupModel() # Setup self model self.drugSource = self.protODE.drugSource self.drugSourceList = self.protODE.drugSourceList self.model = self.protODE.model self.modelList = self.protODE.modelList self.model.deltaT = self.deltaT.get() self.model.setXVar(self.varNameX) self.model.setYVar(self.varNameY) # Get the values to fit x, y = sample.getXYValues(self.varNameX,self.varNameY) print("X= "+str(x)) print("Y= "+str(y)) firstX=x[0] # From [array(...)] to array(...) firstY=y[0] # From [array(...)] to array(...) # Interpret the dose self.protODE.setVarNames(self.varNameX,self.varNameY) self.protODE.setModel(self.model) self.protODE.setTimeRange(sample) sample.interpretDose() self.drugSource.setDoses(sample.parsedDoseList, self.model.t0, self.model.tF) self.protODE.configureSource(self.drugSource) self.model.drugSource = self.drugSource # Prepare the model self.model.setSample(sample) self.calculateParameterUnits(sample) if self.fitting.modelParameterUnits==None: self.fitting.modelParameterUnits = self.parameterUnits # Get the initial parameters if parameterNames==None: parameterNames = self.getParameterNames() parameters0 = [] for parameterName in parameterNames: parameters0.append(float(sample.descriptors[parameterName])) print("Initial solution: %s"%str(parameters0)) print(" ") # Set bounds self.setBounds(sample) # Output object sampleFit = PKPDSampleFitBootstrap() sampleFit.sampleName = sample.sampleName sampleFit.parameters = np.zeros((self.Nbootstrap.get(),len(parameters0)),np.double) sampleFit.xB = [] sampleFit.yB = [] # Bootstrap samples idx = [k for k in range(0,len(firstX))] for n in range(0,self.Nbootstrap.get()): ok = False while not ok: if self.sampleLength.get()>0: lenToUse = self.sampleLength.get() else: lenToUse = len(idx) idxB = sorted(np.random.choice(idx,lenToUse)) xB = [np.asarray([firstX[i] for i in idxB])] yB = [np.asarray([firstY[i] for i in idxB])] print("Bootstrap sample %d"%n) print("X= "+str(xB)) print("Y= "+str(yB)) self.clearXYLists() self.setXYValues(xB, yB) self.parameters = parameters0 optimizer2 = PKPDLSOptimizer(self,fitType) optimizer2.verbose = 0 try: optimizer2.optimize(ftol=1e-4, xtol=1e-4) ok=True except Exception as e: print(e) raise(e) ok=False # Evaluate the quality on the whole data set self.clearXYLists() self.setXYValues(x, y) optimizer2.evaluateQuality() print(optimizer2.optimum) print(" R2 = %f R2Adj=%f AIC=%f AICc=%f BIC=%f"%(optimizer2.R2,optimizer2.R2adj,optimizer2.AIC,\ optimizer2.AICc,optimizer2.BIC)) # Keep this result sampleFit.parameters[n,:] = optimizer2.optimum sampleFit.xB.append(str(xB[0])) sampleFit.yB.append(str(yB[0])) sampleFit.copyFromOptimizer(optimizer2) self.fitting.sampleFits.append(sampleFit) self.fitting.modelParameters = self.getParameterNames() self.fitting.modelDescription = self.getDescription() self.fitting.write(self._getPath("bootstrapPopulation.pkpd"))
[docs] def createOutputStep(self): self._defineOutputs(outputPopulation=self.fitting) self._defineSourceRelation(self.inputODE.get(), self.fitting)
#--------------------------- INFO functions -------------------------------------------- def _summary(self): msg = [] msg.append("Number of bootstrap realizations: %d"%self.Nbootstrap.get()) msg.append("Confidence interval: %f"%self.confidenceInterval.get()) return msg def _validate(self): return []