<pre><code>y_fit = ls.model(t_fit, best_frequency)
</code></pre>
<p>I figured it out. The input t_fit in the code line above is taken as <code>0.0</code> to <code>1.0 DAYS</code> instead of <code>PHASES</code> so the <code>y_fit</code> values are the values of a <code>1.2d</code> orbit but only up to <code>1.0</code> days. If the <code>t_fit</code> range is set instead as (0.0 - 1/best_frequency), we get the full phase of the period orbit. To then plot this full orbit ranging from 0 - 1 over the phase-folded data, the <code>t_fit</code> can be divided by <code>1/best_frequency</code> and it then fits the data correctly.</p>

CC BY-SA 4.0
2021-09-21T12:45:12.547Z
2021-09-21T12:48:28.087Z
2021-09-21T12:48:28.087Z
Mukul Kumar
<p>I am using the astropy Lomb-Scargle periodogram here, <a href="https://docs.astropy.org/en/stable/timeseries/lombscargle.html" rel="nofollow noreferrer">https://docs.astropy.org/en/stable/timeseries/lombscargle.html</a> to get the best period fit using its model() method for simulated data given in their examples. I am then plotting this fit over the phase-folded data using the pyasl.foldAt() method but the period fits are not within the range of 0 - 1 for periods above 1d?</p>
<p>When the periodicity within the simulated data is set to 1d, the best period fit works perfectly as shown below; <a href="https://i.sstatic.net/WXO3El.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/WXO3El.png" alt="Best_period_fit_1d" /></a></p>
<p>But if the period within the simulated data is set to 1.2d, the best period fit is not contained within the 0 - 1 phase range even though the t_fit array given to the ls.model() method is from 0 - 1, plot shown below; <a href="https://i.sstatic.net/oFrCbl.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/oFrCbl.png" alt="Best_period_fit_1.2d" /></a></p>
<p>Below is the code I used.</p>
<pre><code>%matplotlib notebook

import numpy as np
import matplotlib.pyplot as plt
from PyAstronomy import pyasl
from astropy.timeseries import LombScargle

# Quick benchmark of this function on a simple simulated and noisy periodic data!

period = 1.2 # periodicity within the simulated data
rand = np.random.default_rng(42)
t = 100 * rand.random(100)
dy = 0.1 * (1 + rand.random(100))
y = np.sin(2 * np.pi * t * (1/period)) + dy * rand.standard_normal(100)

frequency, power = LombScargle(t, y, dy).autopower()
print('Period at max. power: {}'.format(1/frequency[np.argmax(power)]))

best_frequency = frequency[np.argmax(power)]
t_fit = np.linspace(0.0, 1.0, 1000)
ls = LombScargle(t, y, dy)
y_fit = ls.model(t_fit, best_frequency)

phases = pyasl.foldAt(t, period=1/frequency[np.argmax(power)])

plt.errorbar(phases, y, yerr=dy, fmt='.b')
plt.plot(t_fit, y_fit)
plt.ylabel('Measurement values')
plt.xlabel('Time (d)')
plt.show()
</code></pre>

CC BY-SA 4.0
2021-09-21T11:14:37.193Z
Have you checked whether the frequency corresponding to 1 day is contained in the peridogram? Maybe autopower is failing to include that frequency. You could try to force it to include it by using the options "minimum_frequency" and "maximum_frequency"
Prallax
CC BY-SA 4.0
2021-09-21T11:34:03.050Z
@Prallax, Yes, the periodogram does include the 1d frequency in the periodogram and finds it as its peak, i.e. period with max. power. The issue is that when I use ls.model() to calculate y_fit values for a period of 1.2 days with t_fit ranging from 0-1, it instead just returns the y_fit values actually ranging from 0 - 1.2 days and only returns those up to 1d. If I were to change the t_fit to (0.0, 1.2), it shows the complete phase of 1.2 days but its not fit to the data at all
Mukul Kumar
<exoplanet><astropy>
Incorrect best period fit over phase-folded data
Mukul Kumar
# Incorrect best period fit over phase-folded data

(Asked by: username_0 on 2021-09-21)

I am using the astropy Lomb-Scargle periodogram here, https://docs.astropy.org/en/stable/timeseries/lombscargle.html to get the best period fit using its model() method for simulated data given in their examples. I am then plotting this fit over the phase-folded data using the pyasl.foldAt() method but the period fits are not within the range of 0 - 1 for periods above 1d?
When the periodicity within the simulated data is set to 1d, the best period fit works perfectly as shown below; 
But if the period within the simulated data is set to 1.2d, the best period fit is not contained within the 0 - 1 phase range even though the t_fit array given to the ls.model() method is from 0 - 1, plot shown below; 
Below is the code I used.
%matplotlib notebook

import numpy as np
import matplotlib.pyplot as plt
from PyAstronomy import pyasl
from astropy.timeseries import LombScargle

# Quick benchmark of this function on a simple simulated and noisy periodic data!

period = 1.2 # periodicity within the simulated data
rand = np.random.default_rng(42)
t = 100 * rand.random(100)
dy = 0.1 * (1 + rand.random(100))
y = np.sin(2 * np.pi * t * (1/period)) + dy * rand.standard_normal(100)

frequency, power = LombScargle(t, y, dy).autopower()
print('Period at max. power: {}'.format(1/frequency[np.argmax(power)]))

best_frequency = frequency[np.argmax(power)]
t_fit = np.linspace(0.0, 1.0, 1000)
ls = LombScargle(t, y, dy)
y_fit = ls.model(t_fit, best_frequency)

phases = pyasl.foldAt(t, period=1/frequency[np.argmax(power)])

plt.errorbar(phases, y, yerr=dy, fmt='.b')
plt.plot(t_fit, y_fit)
plt.ylabel('Measurement values')
plt.xlabel('Time (d)')
plt.show()

**Tags:** exoplanet, astropy
> **Comments:**
> **username_1**: Have you checked whether the frequency corresponding to 1 day is contained in the peridogram? Maybe autopower is failing to include that frequency. You could try to force it to include it by using the options "minimum_frequency" and "maximum_frequency"  
> **username_0**: @Prallax, Yes, the periodogram does include the 1d frequency in the periodogram and finds it as its peak, i.e. period with max. power. The issue is that when I use ls.model() to calculate y_fit values for a period of 1.2 days with t_fit ranging from 0-1, it instead just returns the y_fit values actually ranging from 0 - 1.2 days and only returns those up to 1d. If I were to change the t_fit to (0.0, 1.2), it shows the complete phase of 1.2 days but its not fit to the data at all  

---

## 1 answers

### Answer by username_0 on 2021-09-21. Score: 1

y_fit = ls.model(t_fit, best_frequency)

I figured it out. The input t_fit in the code line above is taken as 0.0 to 1.0 DAYS instead of PHASES so the y_fit values are the values of a 1.2d orbit but only up to 1.0 days. If the t_fit range is set instead as (0.0 - 1/best_frequency), we get the full phase of the period orbit. To then plot this full orbit ranging from 0 - 1 over the phase-folded data, the t_fit can be divided by 1/best_frequency and it then fits the data correctly.


---

