diff --git a/bob/learn/em/factor_analysis.py b/bob/learn/em/factor_analysis.py
index 184c5b7c4bfed93ba5471c2e3de0f4a6f0406fb0..1722053250b398edbef8baa88ca4995e9d95968d 100644
--- a/bob/learn/em/factor_analysis.py
+++ b/bob/learn/em/factor_analysis.py
@@ -1136,7 +1136,7 @@ class FactorAnalysisBase(BaseEstimator):
 
         return fn_x.flatten()
 
-    def score_with_array(self, model, data):
+    def score(self, model, data):
         """
         Computes the ISV score using a numpy array as input
 
@@ -1155,7 +1155,7 @@ class FactorAnalysisBase(BaseEstimator):
 
         """
 
-        return self.score(model, self.ubm.transform(data))
+        return self.score_using_stats(model, self.ubm.transform(data))
 
 
 class ISVMachine(FactorAnalysisBase):
@@ -1349,7 +1349,7 @@ class ISVMachine(FactorAnalysisBase):
         """
         return self.enroll([self.ubm.transform(X)], iterations)
 
-    def score(self, latent_z, data):
+    def score_using_stats(self, latent_z, data):
         """
         Computes the ISV score
 
@@ -1839,7 +1839,7 @@ class JFAMachine(FactorAnalysisBase):
 
         return self
 
-    def score(self, model, data):
+    def score_using_stats(self, model, data):
         """
         Computes the JFA score
 
diff --git a/bob/learn/em/test/test_factor_analysis.py b/bob/learn/em/test/test_factor_analysis.py
index fc2ed1595ad98476d80b6df715df6002c7500e44..3e03e9531c4988533ae2a8e580ba9b63460798b4 100644
--- a/bob/learn/em/test/test_factor_analysis.py
+++ b/bob/learn/em/test/test_factor_analysis.py
@@ -392,14 +392,14 @@ def test_JFAMachine():
     model = [y, z]
 
     score_ref = -2.111577181208289
-    score = m.score(model, gs)
+    score = m.score_using_stats(model, gs)
     np.testing.assert_allclose(score, score_ref, atol=eps)
 
     # Scoring with numpy array
     np.random.seed(0)
     X = np.random.normal(loc=0.0, scale=1.0, size=(50, 3))
     score_ref = 2.028009315286946
-    score = m.score_with_array(model, X)
+    score = m.score(model, X)
     np.testing.assert_allclose(score, score_ref, atol=eps)
 
 
@@ -431,7 +431,7 @@ def test_ISVMachine():
 
     # Enrolled model
     latent_z = np.array([3, 4, 1, 2, 0, 1], "float64")
-    score = isv_machine.score(latent_z, gs)
+    score = isv_machine.score_using_stats(latent_z, gs)
     score_ref = -3.280498193082100
     np.testing.assert_allclose(score, score_ref, atol=eps)
 
@@ -439,5 +439,5 @@ def test_ISVMachine():
     np.random.seed(0)
     X = np.random.normal(loc=0.0, scale=1.0, size=(50, 3))
     score_ref = -1.2343813195374242
-    score = isv_machine.score_with_array(latent_z, X)
+    score = isv_machine.score(latent_z, X)
     np.testing.assert_allclose(score, score_ref, atol=eps)