pcVPC with Stratification

1

Hi,

Currently I can see that

ArgumentError: Prediction correction is not currently supported when stratifying the VPC

is there a workaround, were I can convert the vpc object to a dataframe and then use it for plotting.

Thanks for the help.

Best,

Parssh

2

Dear Paarsh - Generally, it is not recommended to stratify when you perform prediction correction as it defeats the purpose of prediction correction.

To answer your other question on how to download the data from a vpc run: any vpc object that you create has fields within it which store both the simulated data and the observed data. You can easily access these fields by calling them e.g. myvpcobject.data

We are working on newer and easier ways to extract the raw data out of a vpc object which we will post as a follow up to here.

Best,
Vijay

3

Hi Vijay,

Thanks for the prompt reply. I totally understand that the PRED calculation will account for the covariate effect as it will be calculated accordingly. But this is a requirement from the regulatory agency. I was able to extract the quantiles using the code below. But the dataframe does not contain the covariates, is there a way to get them as well.

DataFrame(pk_vpc_cp_fit17.simulated_quantiles)

Best,

Parssh

4

Dear Parssh,

The covariates are not part of the output because the smoothing is done across the chosen “x-axis” through the covariates keyword argument to vpc. If you were able to stratify, the strata would be part of the data and you could plot it, but as the error message states, it is not possible.

The reason why we present you with the error message is that in a previous version of Pumas it was allowed, but the output was incorrect. We have significant improvements coming to the VPC workflow and then we will re-enable this combination of inputs (prediction correction + stratification).

That said, I think we can find a solution. Can you describe the covariate you with to stratify on? What is its names and what is the nature of the unique values? If you cannot disclose the exact information that is fine, but is it a binary, a limited number of categorical values, or ?

Best,
Patrick

5

Thanks Patrick. So the covariates are X= Two levels [1,2] and Y= Two levels [1,2]. I need the stratification to have a combination of [1,1], [1,2] & [2,2]. The x axis will be :tad.

Best,

Parssh

6

Okay. You will need to subset your population. You can either create these subpopulations already from the data frames and with read_pumas, or we can filter later. The vpc function has a function that accepts a model, population, and coefficients. This following has covariates so illustrates the idea, but it is not a continuous outcome so you cannot use prediction correction here :slight_smile:

using PumasPlots
using PumasPlots.Pumas
using Pumas.CSV
using Pumas.DataFrames
using Random
using StableRNGs
import PharmaDatasets

df = PharmaDatasets.dataset("pumas/pain_remed")
pop = read_pumas(
    df;
    observations = [:dv],
    covariates = [:arm, :dose, :conc, :painord, :remed],
    event_data = false,
)
model = @model begin
    @metadata begin
        desc = "Logistic Regression Model"
    end
    @param begin
        "Intercept"
        intercept ∈ RealDomain(init = 0.001)
        "Slope"
        tvslope ∈ RealDomain(init = 0.0001)
        """
         - Ω
        """
        Ω ∈ VectorDomain(1)
    end
    @random begin
        η ~ MvNormal(Ω)
    end
    @covariates begin
        "Study Arm"
        arm
        "Dose (mg)"
        dose
    end
    @pre begin
        rx = dose > 0 ? 1 : 0
        slope = tvslope * rx
        logit = intercept + slope + η[1]
    end
    @derived begin
        "Pain Relief"
        dv ~ @. Bernoulli(logistic(logit))
    end
end
param = (intercept = 0.001, tvslope = 0.0001, Ω = [1.0])
fit1 = fit(model, pop, param, Pumas.FOCE(), ensemblealg = EnsembleThreads())
vpc1 = vpc(
    model,
    pop,
    coef(fit1),
    samples = 200,
    stratify_by = [:arm, :dose],
    ensemblealg = EnsembleThreads(),
    rng = StableRNG(1),
)
vpc_plot(vpc1)

pop_arm1_dose5 = filter(x -> x.covariates(0.0).arm == 1 && x.covariates(0.0).dose == 5, pop)
#julia> pop_arm1[1].covariates(0)
#(arm = 1,
# dose = 5,
# conc = 0.0,
# painord = 2,
# remed = 0,)

vpc_arm1_dose5 = vpc(
    model,
    pop_arm1_dose5,
    coef(fit1),
    samples = 200,
    stratify_by = [:arm, :dose],
    ensemblealg = EnsembleThreads(),
    rng = StableRNG(1),
#    prediction_correction=true,
)

vpc_plot(vpc_arm1_dose5)
7

Thanks a lot Patrick. That was helpful.

8

If you need further help, let me know!

1 Like