Repository: hgrecco/pint-pandas

Issue:
Very long runtimes for  DataFrames with units
Using pint dtypes increases run time significantly compared to standard pandas DataFrames.

Example:

```python

import numpy as np
import pandas as pd
import pint_pandas


def make_df(size, pint_units=True, dtype=float):
    if pint_units:
        dist_unit ='pint[m]'
        time_unit = 'pint[s]'
    else:
        dist_unit = dtype
        time_unit = dtype
    return pd.DataFrame(
        {'distance': pd.Series(np.arange(1, size + 1, dtype=dtype), dtype=dist_unit),
         'time':  pd.Series(np.arange(1, size + 1, dtype=dtype), dtype=time_unit)
        }
    )


n = 1_000_000
df_pint = make_df(n)
df = make_df(n, pint_units=False)
```

Now, time in a Notebook:

```python
tp = %timeit -o df_pint['speed'] = df_pint['distance'] / df_pint[('time')]
822 ms ± 12 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

t = %timeit -o df['speed'] = df.distance * df.time
3.61 ms ± 153 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

tp.average / t.average
227.46385324897062
```

The unit version takes more than 200x as long as the DataFrame without units. Am I doing something wrong? The computation should be done by NumPy in the background. Therefore, I would expect very litte performance differences. 

Discussion:

Patch:
diff --git a/CHANGES b/CHANGES
index 12c2d3b..746db0f 100644
--- a/CHANGES
+++ b/CHANGES
@@ -4,7 +4,9 @@ pint-pandas Changelog
 0.3 (unreleased)
 ----------------
 
-- Nothing changed yet.
+- Fixed a bug which caused the creation of PintArray (and even more so all binary mathematical operations)
+  to convert input arrays to lists of Quantity, thereby causing a huge slowdown. #80
+- Fixed a bug where listlike operands in binary operations were not correctly converted to Quantity arrays. #86
 
 
 0.2 (2021-03-23)
diff --git a/pint_pandas/pint_array.py b/pint_pandas/pint_array.py
index 35f2fad..50505dd 100644
--- a/pint_pandas/pint_array.py
+++ b/pint_pandas/pint_array.py
@@ -15,6 +15,7 @@ from pandas.api.extensions import (
 )
 from pandas.api.types import is_integer, is_list_like, is_object_dtype, is_string_dtype
 from pandas.compat import set_function_name
+from pandas.core import nanops
 from pandas.core.arrays.base import ExtensionOpsMixin
 from pandas.core.indexers import check_array_indexer
 from pint import compat, errors
@@ -189,13 +190,19 @@ class PintArray(ExtensionArray, ExtensionOpsMixin):
         if not isinstance(dtype, PintType):
             dtype = PintType(dtype)
         self._dtype = dtype
-        if len(values) > 0 and all([not isinstance(x, float) for x in values]):
-            data_dtype = next(x for x in values if not isinstance(x, float))
+        if not isinstance(values, np.ndarray):
+            values = np.array(values, copy=copy)
+            copy = False
+        if not np.issubdtype(values.dtype, np.floating):
             warnings.warn(
-                f"pint-pandas does not support magnitudes of {type(data_dtype)}. Converting magnitudes to float.",
+                f"pint-pandas does not support magnitudes of {values.dtype}. Converting magnitudes to float.",
                 category=RuntimeWarning,
             )
-        self._data = np.array(values, float, copy=copy)
+            values = values.astype(float)
+            copy = False
+        if copy:
+            values = values.copy()
+        self._data = values
         self._Q = self.dtype.ureg.Quantity
 
     @property
@@ -565,38 +572,31 @@ class PintArray(ExtensionArray, ExtensionOpsMixin):
 
         def _binop(self, other):
             def validate_length(obj1, obj2):
-                # validates length and converts to listlike
+                # validates length
+                # CHANGED: do not convert to listlike (why should we? pint.Quantity is perfecty able to handle that...)
                 try:
-                    if len(obj1) == len(obj2):
-                        return obj2
-                    else:
+                    if len(obj1) != len(obj2):
                         raise ValueError("Lengths must match")
                 except TypeError:
-                    return [obj2] * len(obj1)
+                    pass
 
             def convert_values(param):
                 # convert to a quantity or listlike
                 if isinstance(param, cls):
                     return param.quantity
-                elif isinstance(param, _Quantity):
+                elif isinstance(param, (_Quantity, _Unit)):
                     return param
                 elif is_list_like(param) and isinstance(param[0], _Quantity):
-                    return type(param[0])([p.magnitude for p in param], param[0].units)
+                    units = param[0].units
+                    return type(param[0])([p.m_as(units) for p in param], units)
                 else:
                     return param
 
             if isinstance(other, (Series, DataFrame)):
                 return NotImplemented
             lvalues = self.quantity
-            other = validate_length(lvalues, other)
+            validate_length(lvalues, other)
             rvalues = convert_values(other)
-            # Pint quantities may only be exponented by single values, not arrays.
-            # Reduce single value arrays to single value to allow power ops
-            if isinstance(rvalues, _Quantity):
-                if len(set(np.array(rvalues.data))) == 1:
-                    rvalues = rvalues[0]
-            elif len(set(np.array(rvalues))) == 1:
-                rvalues = rvalues[0]
             # If the operator is not defined for the underlying objects,
             # a TypeError should be raised
             res = op(lvalues, rvalues)
@@ -699,7 +699,7 @@ class PintArray(ExtensionArray, ExtensionOpsMixin):
             value = [item.to(self.units).magnitude for item in value]
         return arr.searchsorted(value, side=side, sorter=sorter)
 
-    def _reduce(self, name, skipna=True, **kwds):
+    def _reduce(self, name, **kwds):
         """
         Return a scalar result of performing the reduction operation.
 
@@ -723,24 +723,30 @@ class PintArray(ExtensionArray, ExtensionOpsMixin):
         ------
         TypeError : subclass does not define reductions
         """
+
         functions = {
-            "all": all,
-            "any": any,
-            "min": min,
-            "max": max,
-            "sum": sum,
-            "mean": np.mean,
-            "median": np.median,
+            "any": nanops.nanany,
+            "all": nanops.nanall,
+            "min": nanops.nanmin,
+            "max": nanops.nanmax,
+            "sum": nanops.nansum,
+            "mean": nanops.nanmean,
+            "median": nanops.nanmedian,
+            "std": nanops.nanstd,
+            "var": nanops.nanvar,
+            "sem": nanops.nansem,
+            "kurt": nanops.nankurt,
+            "skew": nanops.nanskew,
         }
         if name not in functions:
             raise TypeError(f"cannot perform {name} with type {self.dtype}")
 
-        if skipna:
-            quantity = self.dropna().quantity
-        else:
-            quantity = self.quantity
-
-        return functions[name](quantity)
+        result = functions[name](self._data, **kwds)
+        if name in {"all", "any", "kurt", "skew"}:
+            return result
+        if name == "var":
+            return self._Q(result, self.units ** 2)
+        return self._Q(result, self.units)
 
 
 PintArray._add_arithmetic_ops()


Test patch:
diff --git a/pint_pandas/testsuite/test_issues.py b/pint_pandas/testsuite/test_issues.py
new file mode 100644
index 0000000..67b4ff4
--- /dev/null
+++ b/pint_pandas/testsuite/test_issues.py
@@ -0,0 +1,73 @@
+import time
+
+import numpy as np
+import pandas as pd
+import pytest
+
+import pint_pandas as ppi
+from pint_pandas import PintArray
+
+ureg = ppi.PintType.ureg
+
+
+class TestIssue80:
+    @staticmethod
+    def _timeit(fun, n_runs=5):
+        run_time = []
+        for k in range(n_runs):
+            t_start = time.monotonic_ns()
+            fun()
+            t_end = time.monotonic_ns()
+            run_time.append(t_end - t_start)
+        return np.median(run_time) * ureg.ns
+
+    @staticmethod
+    def _make_df(size, pint_units=True, dtype=float):
+        if pint_units:
+            dist_unit = "pint[m]"
+            time_unit = "pint[s]"
+        else:
+            dist_unit = dtype
+            time_unit = dtype
+        return pd.DataFrame(
+            {
+                "distance": pd.Series(
+                    np.arange(1, size + 1, dtype=dtype), dtype=dist_unit
+                ),
+                "time": pd.Series(np.arange(1, size + 1, dtype=dtype), dtype=time_unit),
+            }
+        )
+
+    def test_div(self):
+        n = 1_000_000
+        df_pint = self._make_df(n)
+        df = self._make_df(n, pint_units=False)
+
+        tp = self._timeit(lambda: df_pint["distance"] / df_pint["time"]).to("ms")
+        t = self._timeit(lambda: df["distance"] / df["time"]).to("ms")
+
+        assert tp < 5 * t
+
+    @pytest.mark.parametrize(
+        "reduction",
+        ["min", "max", "sum", "mean", "median"],
+    )
+    def test_reductions(self, reduction):
+        # before the fix, those reductions could be *very* slow. Fail early.
+        for n in [10_000, 1_000_000]:
+            s_pint = self._make_df(n)["time"]
+            s = self._make_df(n, pint_units=False)["time"]
+
+            tp = self._timeit(getattr(s_pint, reduction)).to("ms")
+            t = self._timeit(getattr(s, reduction)).to("ms")
+
+            assert tp < 5 * t
+
+
+def test_issue_86():
+    a = PintArray([1, 2], ureg.m)
+    b_listlike = [1 * ureg.km, 1 * ureg.m]
+    units = b_listlike[0].units
+    b_pa = PintArray([v.m_as(units) for v in b_listlike], units)
+
+    assert np.all(a + b_listlike == a + b_pa)
diff --git a/pint_pandas/testsuite/test_pandas_interface.py b/pint_pandas/testsuite/test_pandas_interface.py
index a3f4f64..8b5093a 100644
--- a/pint_pandas/testsuite/test_pandas_interface.py
+++ b/pint_pandas/testsuite/test_pandas_interface.py
@@ -1,5 +1,6 @@
 import itertools
 import operator
+import warnings
 from os.path import dirname, join
 
 import numpy as np
@@ -170,11 +171,12 @@ _all_numeric_reductions = [
     "min",
     "mean",
     # "prod",
-    # "std",
-    # "var",
+    "std",
+    "var",
     "median",
-    # "kurt",
-    # "skew",
+    "sem",
+    "kurt",
+    "skew",
 ]
 
 
@@ -622,10 +624,44 @@ class TestNumericReduce(base.BaseNumericReduceTests):
         expected_m = getattr(pd.Series(s.values.quantity._magnitude), op_name)(
             skipna=skipna
         )
-        expected_u = s.values.quantity.units
-        expected = ureg.Quantity(expected_m, expected_u)
+        if op_name in {"kurt", "skew"}:
+            expected_u = None
+        elif op_name in {"var"}:
+            expected_u = s.values.quantity.units ** 2
+        else:
+            expected_u = s.values.quantity.units
+        if expected_u is not None:
+            expected = ureg.Quantity(expected_m, expected_u)
+        else:
+            expected = expected_m
         assert result == expected
 
+    @pytest.mark.parametrize("skipna", [True, False])
+    def test_reduce_scaling(self, data, all_numeric_reductions, skipna):
+        """Make sure that the reductions give the same physical result independent of the unit representation.
+
+        This verifies that the result units are sensible.
+        """
+        op_name = all_numeric_reductions
+        s_nm = pd.Series(data)
+        # Attention: `mm` is fine here, but with `m`, the magnitudes become so small
+        # that pandas discards them in the kurtosis calculation, leading to different results.
+        s_mm = pd.Series(PintArray.from_1darray_quantity(data.quantity.to(ureg.mm)))
+
+        # min/max with empty produce numpy warnings
+        with warnings.catch_warnings():
+            warnings.simplefilter("ignore", RuntimeWarning)
+            r_nm = getattr(s_nm, op_name)(skipna=skipna)
+            r_mm = getattr(s_mm, op_name)(skipna=skipna)
+            if isinstance(r_nm, ureg.Quantity):
+                # convert both results to the same units, then take the magnitude
+                v_nm = r_nm.m_as(r_mm.units)
+                v_mm = r_mm.m
+            else:
+                v_nm = r_nm
+                v_mm = r_mm
+            assert np.isclose(v_nm, v_mm, rtol=1e-3), f"{r_nm} == {r_mm}"
+
 
 class TestBooleanReduce(base.BaseBooleanReduceTests):
     def check_reduce(self, s, op_name, skipna):
@@ -989,6 +1025,11 @@ comparative_ops = [
     operator.gt,
 ]
 
+unit_ops = [
+    operator.mul,
+    operator.truediv,
+]
+
 
 class TestPintArrayQuantity(QuantityTestCase):
     FORCE_NDARRAY = True
@@ -1042,12 +1083,18 @@ class TestPintArrayQuantity(QuantityTestCase):
             ureg.Quantity([7.0, np.nan]),
         ]
 
+        us = [ureg.m]
+
         for a_pint, a_pint_array in zip(a_pints, a_pint_arrays):
             for b in bs:
                 for op in arithmetic_ops:
                     test_op(a_pint, a_pint_array, b)
                 for op in comparative_ops:
                     test_op(a_pint, a_pint_array, b, coerce=False)
+            # also test for operations involving units
+            for b in us:
+                for op in unit_ops:
+                    test_op(a_pint, a_pint_array, b)
 
     def test_mismatched_dimensions(self):
         x_and_ys = [