# pcVPC with Stratification

**URL:** <https://discourse.pumas.ai/t/pcvpc-with-stratification/1133>\
**Category:** Uncategorized\
**Created:** [February 9, 2025, 6:05pm UTC](https://discourse.pumas.ai/t/pcvpc-with-stratification/1133 "2025-02-09T18:05:22Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![parssh.mehta](https://avatars.discourse-cdn.com/v4/letter/p/8c91f0/32.png) [@parssh.mehta](https://discourse.pumas.ai/u/parssh.mehta)\
**Post date:** [February 9, 2025, 6:05pm UTC](https://discourse.pumas.ai/t/pcvpc-with-stratification/1133/1 "2025-02-09T18:05:22Z")

</div>

Hi,

Currently I can see that

```auto
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

---

<div class="post-metadata">

**Author:** ![vijay](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/vijay/32/96_2.png) [@vijay](https://discourse.pumas.ai/u/vijay)\
**Post date:** [February 10, 2025, 12:18am UTC](https://discourse.pumas.ai/t/pcvpc-with-stratification/1133/2 "2025-02-10T00:18:28Z")

</div>

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

---

<div class="post-metadata">

**Author:** ![parssh.mehta](https://avatars.discourse-cdn.com/v4/letter/p/8c91f0/32.png) [@parssh.mehta](https://discourse.pumas.ai/u/parssh.mehta)\
**Post date:** [February 10, 2025, 3:08pm UTC](https://discourse.pumas.ai/t/pcvpc-with-stratification/1133/3 "2025-02-10T15:08:50Z")

</div>

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.

```auto
DataFrame(pk_vpc_cp_fit17.simulated_quantiles)

```

Best,

Parssh

---

<div class="post-metadata">

**Author:** ![patrick](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/patrick/32/206_2.png) [@patrick](https://discourse.pumas.ai/u/patrick)\
**Post date:** [February 11, 2025, 8:02am UTC](https://discourse.pumas.ai/t/pcvpc-with-stratification/1133/4 "2025-02-11T08:02:05Z")

</div>

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

---

<div class="post-metadata">

**Author:** ![parssh.mehta](https://avatars.discourse-cdn.com/v4/letter/p/8c91f0/32.png) [@parssh.mehta](https://discourse.pumas.ai/u/parssh.mehta)\
**Post date:** [February 11, 2025, 2:19pm UTC](https://discourse.pumas.ai/t/pcvpc-with-stratification/1133/5 "2025-02-11T14:19:27Z")

</div>

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

---

<div class="post-metadata">

**Author:** ![patrick](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/patrick/32/206_2.png) [@patrick](https://discourse.pumas.ai/u/patrick)\
**Post date:** [February 14, 2025, 3:12pm UTC](https://discourse.pumas.ai/t/pcvpc-with-stratification/1133/6 "2025-02-14T15:12:51Z")

</div>

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 🙂

```auto
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)

```

---

<div class="post-metadata">

**Author:** ![parssh.mehta](https://avatars.discourse-cdn.com/v4/letter/p/8c91f0/32.png) [@parssh.mehta](https://discourse.pumas.ai/u/parssh.mehta)\
**Post date:** [February 14, 2025, 11:57pm UTC](https://discourse.pumas.ai/t/pcvpc-with-stratification/1133/7 "2025-02-14T23:57:39Z")

</div>

Thanks a lot Patrick. That was helpful.

---

<div class="post-metadata">

**Author:** ![patrick](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/patrick/32/206_2.png) [@patrick](https://discourse.pumas.ai/u/patrick)\
**Post date:** [February 18, 2025, 10:38am UTC](https://discourse.pumas.ai/t/pcvpc-with-stratification/1133/8 "2025-02-18T10:38:26Z")

</div>

If you need further help, let me know!