diff --git a/quapy/data/base.py b/quapy/data/base.py
index 6db182a..9bdf135 100644
--- a/quapy/data/base.py
+++ b/quapy/data/base.py
@@ -44,7 +44,7 @@ class LabelledCollection:
 
     @property
     def index(self):
-        if self._index is None:
+        if not hasattr(self, '_index') or self._index is None:
             self._index = {class_: np.arange(len(self))[self.labels == class_] for class_ in self.classes_}
         return self._index
 
diff --git a/quapy/method/confidence.py b/quapy/method/confidence.py
index dd9e05f..7423545 100644
--- a/quapy/method/confidence.py
+++ b/quapy/method/confidence.py
@@ -13,6 +13,7 @@ from abc import ABC, abstractmethod
 from scipy.special import softmax, factorial
 import copy
 from functools import lru_cache
+from pathlib import Path
 
 """
 This module provides implementation of different types of confidence regions, and the implementation of Bootstrap
@@ -613,7 +614,7 @@ class PQ(AggregativeSoftQuantifier, BinaryAggregativeQuantifier):
         if num_samples <= 0:
             raise ValueError(f'parameter {num_samples=} must be a positive integer')
 
-        if _bayesian.DEPENDENCIES_INSTALLED is False:
+        if not _bayesian.DEPENDENCIES_INSTALLED:
             raise ImportError("Auxiliary dependencies are required. "
                               "Run `$ pip install quapy[bayes]` to install them.")
 
@@ -624,7 +625,9 @@ class PQ(AggregativeSoftQuantifier, BinaryAggregativeQuantifier):
         self.num_samples = num_samples
         self.region = region
         self.stan_seed = stan_seed
-        with open('quapy/method/stan/pq.stan', 'r') as f:
+        # with open('quapy/method/stan/pq.stan', 'r') as f:
+        stan_path = Path(__file__).resolve().parent / "stan" / "pq.stan"
+        with stan_path.open("r") as f:
             self.stan_code = str(f.read())
 
     def aggregation_fit(self, classif_predictions, labels):
@@ -649,24 +652,17 @@ class PQ(AggregativeSoftQuantifier, BinaryAggregativeQuantifier):
         self.pos_hist = np.bincount(bin_indices[pos_mask], minlength=self.n_bins)
         self.neg_hist = np.bincount(bin_indices[neg_mask], minlength=self.n_bins)
 
-
     def aggregate(self, classif_predictions):
         Px_test = classif_predictions[:, self.pos_label]
         test_hist, _ = np.histogram(Px_test, bins=self.bin_limits)
-        self.prev_distribution = _bayesian.pq_stan(self.stan_code, 
-                                                            self.n_bins, 
-                                                            self.pos_hist, 
-                                                            self.neg_hist, 
-                                                            test_hist, 
-                                                            self.num_samples, 
-                                                            self.num_warmup, 
-                                                            self.stan_seed)
+        self.prev_distribution = _bayesian.pq_stan(
+            self.stan_code, self.n_bins, self.pos_hist, self.neg_hist, test_hist,
+            self.num_samples, self.num_warmup, self.stan_seed
+        )
         return F.as_binary_prevalence(self.prev_distribution.mean())
-        
 
     def predict_conf(self, instances, confidence_level=None) -> (np.ndarray, ConfidenceRegionABC):
-        classif_predictions = self.classify(instances)
-        point_estimate = self.aggregate(classif_predictions)
+        point_estimate = self.predict(instances)
         samples = self.prev_distribution
         region = WithConfidenceABC.construct_region(samples, confidence_level=confidence_level, method=self.region)
         return point_estimate, region
diff --git a/quapy/method/stan/pq.stan b/quapy/method/stan/pq.stan
index b759ddf..a5bc52a 100644
--- a/quapy/method/stan/pq.stan
+++ b/quapy/method/stan/pq.stan
@@ -35,4 +35,5 @@ model {
 generated quantities {
   real<lower=0,upper=1> prev;
   prev = sum( binomial_rng(test, 1 / ( 1 + (p_neg./p_pos) *(1-prev_prior)/prev_prior ) ) ) / n_test;
-}
\ No newline at end of file
+}
+
diff --git a/setup.py b/setup.py
index ba5f205..c058f9d 100644
--- a/setup.py
+++ b/setup.py
@@ -124,7 +124,7 @@ setup(
     # Similar to `install_requires` above, these must be valid existing
     # projects.
     extras_require={  # Optional
-       'bayes': ['jax', 'jaxlib', 'numpyro'],
+       'bayes': ['jax', 'jaxlib', 'numpyro', 'pystan'],
        'neural': ['torch'],
        'tests': ['certifi'],
        'docs' : ['sphinx-rtd-theme', 'myst-parser'],