FLRW.luminosity_distance returns unexpected type

[mj-will]

Description

When calling FLRM.luminosity_distance with an instance of pandas.Series as the input, the output is another pandas.Series rather than an instance of astropy.units.Quantity. This is inconsistent with other similar methods, such as comoving_distance or comoving_volume.

Expected behavior

FLRM.luminosity_distance should return an instance of astropy.units.Quantity.

How to Reproduce

1. Install astropy, scipy and pandas using conda
2. Run the code snippet below and check the types of the outputs

from astropy.cosmology import WMAP9 as cosmo
import pandas as pd
z = pd.Series([1, 2])
dl = cosmo.luminosity_distance(z)
print(type(dl))
dc = cosmo.comoving_distance(z)
print(type(dc))
vc = cosmo.comoving_volume(z)
print(type(vc))

On my system, this outputs

<class 'pandas.core.series.Series'>
<class 'astropy.units.quantity.Quantity'>
<class 'astropy.units.quantity.Quantity'>

Versions

Linux-3.10.0-1160.71.1.el7.x86_64-x86_64-with-glibc2.17
Python 3.10.13 (main, Sep 11 2023, 13:44:35) [GCC 11.2.0]
astropy 5.3.4
Numpy 1.26.0
pyerfa 2.0.0.1
Scipy 1.11.3
# Matplotlib is not installed

[github-actions]

Welcome to Astropy 👋 and thank you for your first issue!

A project member will respond to you as soon as possible; in the meantime, please double-check the guidelines for submitting issues and make sure you've provided the requested details.

GitHub issues in the Astropy repository are used to track bug reports and feature requests; If your issue poses a question about how to use Astropy, please instead raise your question in the Astropy Discourse user forum and close this issue.

If you feel that this issue has not been responded to in a timely manner, please send a message directly to the development mailing list. If the issue is urgent or sensitive in nature (e.g., a security vulnerability) please send an e-mail directly to the private e-mail feedback@astropy.org.

[pllim]

@nstarman , is this behavior intentional? Please advise. Thanks!

[nstarman]

This is an issue. But it actually an issue for units. @mhvk might want to see this as well.
The problem is that a Pandas Series * Quantity does not return a Quantity (and a Quantity * Series = Series)!

>>> import pandas as pd
>>> import astropy.units as u

>>> z = pd.Series([1, 2])
>>> d = u.Quantity([10, 20], u.Mpc)

>>> z * d
0    10.0
1    40.0
dtype: float64

[nstarman]

I found #11247, so this is a known problem.
We could change aszarr

astropy/astropy/cosmology/_utils.py

Lines 72 to 75 in e42dfc0

	elif hasattr(z, "shape"): # ducktypes NumPy array
	if hasattr(z, "unit"): # Quantity Column
	return (z << cu.redshift).value # for speed only use enabled equivs
	return z

to check for pandas and convert to numpy?

[nstarman]

But the more optimal solution would be to have Quantity correctly handle interactions with Pandas.

[mhvk]

See my comment in the earlier issue: #11247 (comment)
We cannot change the fact that Series * unit does not work because our code is never called -- that really needs a change in pandas... Given that that very basic way of converting Series into a Quantity cannot be made to work, I don't see much point in changing Quantity at all -- as is, the initializer works properly.

Obviously, one can make changes to specific methods in cosmology if you wish - in a way, the problem here is that in trying to get better performance, one skips making Quantity and therefore passes in Series where they do not belong. You could also add something like np.asanyarray and be done with it, but then of course you also prevent well-behaved duck-types from working (if any exist, not sure! maybe dask array these days?).

Anyway, my bottom-line suggestion would be to raise an issue with pandas.

[nstarman]

We cannot change the fact that Series * unit does not work because our code is never called

The issue here is not Series * unit, which I agree is problematic but we cannot address, but Quantity * Series, which we can!

class Quantity:

    def __mul__(self, other):
        if other.__module__.startswith("pandas"):  # no pandas dependency
            return self * np.array(other)  # now we get a Quantity

Pandas is sufficiently common I can see it being worthwhile special-casing, but if we don't want to do that in Quantity, I can handle it in aszarr here.

[mhvk]

I guess we can do that in principle deeper inside, in __array_ufunc__, perhaps by looking for an __array__ method. But it should be in a try/except, with the check in the except part - there shouldn't be any performance regressions for this!

[pllim]

Rather than patching Quantity, which as you stated, is a separate issue and requires upstream fix, is it possible to apply a more hacky but local patch just for that cosmology API? Then in the other larger issue, we can add a reminder to undo the cosmology one when/if we ever get to it.

[nstarman]

Yeah. In aszarr which processes each input we can do essentially the same thing I suggested for Quantity. I can push a PR soon.