Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

LombScargle.false_alarm_probability fails for Quantities #9236

Closed
naamach opened this issue Sep 15, 2019 · 3 comments · Fixed by #9246
Closed

LombScargle.false_alarm_probability fails for Quantities #9236

naamach opened this issue Sep 15, 2019 · 3 comments · Fixed by #9246

Comments

@naamach
Copy link

naamach commented Sep 15, 2019

Trying to calculate the Lomb Scragle false alarm probability while using Quantities with units throws an error while executing fap_baluev, since the exponent argument has units.
(It works ok when using method='naive')

Using Ubuntu 18.04.3, python 3.7.4, astropy 3.1.2, numpy 1.16.2.

Minimal working example:

import numpy as np
import astropy.units as u
from astropy.timeseries import LombScargle

t = 100 * np.random.rand(10) * u.day
dy = 0.1 * u.mag
y = (np.sin(2 * np.pi * t.value) + dy.value * np.random.randn(10)) * u.mag

ls = LombScargle(t, y, dy)
frequency, power = ls.autopower()
ls.false_alarm_probability(power.max())

EDIT: Syntax highlighting

@bsipocz
Copy link
Member

bsipocz commented Sep 16, 2019

Here is the error in question:

---------------------------------------------------------------------------
UnitConversionError                       Traceback (most recent call last)
~/munka/devel/astropy/astropy/units/quantity_helper/helpers.py in get_converter(from_unit, to_unit)
     31     try:
---> 32         scale = from_unit._to(to_unit)
     33     except UnitsError:

~/munka/devel/astropy/astropy/units/core.py in _to(self, other)
    952         raise UnitConversionError(
--> 953             f"'{self!r}' is not a scaled version of '{other!r}'")
    954 

UnitConversionError: 'Unit("1 / d")' is not a scaled version of 'Unit(dimensionless)'

During handling of the above exception, another exception occurred:

UnitConversionError                       Traceback (most recent call last)
~/munka/devel/astropy/astropy/units/quantity_helper/helpers.py in helper_dimensionless_to_dimensionless(f, unit)
    146     try:
--> 147         return ([get_converter(unit, dimensionless_unscaled)],
    148                 dimensionless_unscaled)

~/munka/devel/astropy/astropy/units/quantity_helper/helpers.py in get_converter(from_unit, to_unit)
     34         return from_unit._apply_equivalencies(
---> 35                 from_unit, to_unit, get_current_unit_registry().equivalencies)
     36     except AttributeError:

~/munka/devel/astropy/astropy/units/core.py in _apply_equivalencies(self, unit, other, equivalencies)
    889             "{} and {} are not convertible".format(
--> 890                 unit_str, other_str))
    891 

UnitConversionError: '1 / d' (frequency) and '' (dimensionless) are not convertible

During handling of the above exception, another exception occurred:

UnitTypeError                             Traceback (most recent call last)
<ipython-input-1-4f85c47ebb3e> in <module>
      9 ls = LombScargle(t, y, dy)
     10 frequency, power = ls.autopower()
---> 11 ls.false_alarm_probability(power.max())

~/munka/devel/astropy/astropy/timeseries/periodograms/lombscargle/core.py in false_alarm_probability(self, power, method, samples_per_peak, nyquist_factor, minimum_frequency, maximum_frequency, method_kwds)
    625                                                    normalization=self.normalization,
    626                                                    method=method,
--> 627                                                    method_kwds=method_kwds)
    628 
    629     def false_alarm_level(self, false_alarm_probability, method='baluev',

~/munka/devel/astropy/astropy/timeseries/periodograms/lombscargle/_statistics.py in false_alarm_probability(Z, fmax, t, y, dy, normalization, method, method_kwds)
    433     method_kwds = method_kwds or {}
    434 
--> 435     return method(Z, fmax, t, y, dy, normalization, **method_kwds)
    436 
    437 

~/munka/devel/astropy/astropy/timeseries/periodograms/lombscargle/_statistics.py in fap_baluev(Z, fmax, t, y, dy, normalization)
    327     # result is 1 - (1 - fap_s) * np.exp(-tau)
    328     # this is much more precise for small numbers
--> 329     return -np.expm1(-tau) + fap_s * np.exp(-tau)
    330 
    331 

~/munka/devel/astropy/astropy/units/quantity.py in __array_ufunc__(self, function, method, *inputs, **kwargs)
    445         # consistent units between two inputs (e.g., in np.add) --
    446         # and the unit of the result (or tuple of units for nout > 1).
--> 447         converters, unit = converters_and_unit(function, method, *inputs)
    448 
    449         out = kwargs.get('out', None)

~/munka/devel/astropy/astropy/units/quantity_helper/converters.py in converters_and_unit(function, method, *args)
    164 
    165         # Determine possible conversion functions, and the result unit.
--> 166         converters, result_unit = ufunc_helper(function, *units)
    167 
    168         if any(converter is False for converter in converters):

~/munka/devel/astropy/astropy/units/quantity_helper/helpers.py in helper_dimensionless_to_dimensionless(f, unit)
    150         raise UnitTypeError("Can only apply '{}' function to "
    151                             "dimensionless quantities"
--> 152                             .format(f.__name__))
    153 
    154 

UnitTypeError: Can only apply 'expm1' function to dimensionless quantities

@pllim
Copy link
Member

pllim commented Sep 16, 2019

tau that is being passed into np.expm1 is:

2.0394872747734722 1 / d

where d is day.

@bsipocz bsipocz added the Bug label Sep 16, 2019
@bsipocz
Copy link
Member

bsipocz commented Sep 16, 2019

note: this only affects the baluev and davies methods, both bootstrap and naive works, as well as the single (that is not documented).

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging a pull request may close this issue.

3 participants