Since mostly I play with financial data and data science concepts, I think it’s a good idea to create a series called Data science in investing in which I share the small things I learn.

For today we have Data science in investing 102 at this link, enjoy: https://medium.com/@sorin.tudor/data-science-in-investing-102-54dd5b0550e3

This is not my script, but I needed to share this because it’s impressive.

Here is the generated script by Gemini to check if the index is overvalued

import requests
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import os

# ==========================================
# 1. CONFIGURATION
# ==========================================
API_KEY = os.environ.get('FMP_KEY')  # <--- REPLACE THIS
BASE_URL = 'https://financialmodelingprep.com/api/v3'

# Top STOXX 50 Heavyweights (approx 60-70% of index weight)
# We use a representative subset to keep the script fast.
# For production, use the full 50 list.
stoxx_tickers = [
    # --- FRANCE (approx. 18) ---
    "MC.PA",     # LVMH Moet Hennessy Louis Vuitton
    "OR.PA",     # L'Oreal
    "RMS.PA",    # Hermes International
    "TTE.PA",    # TotalEnergies
    "SAN.PA",    # Sanofi (Note: Do not confuse with SAN.MC)
    "AIR.PA",    # Airbus
    "SU.PA",     # Schneider Electric
    "AI.PA",     # Air Liquide
    "BNP.PA",    # BNP Paribas
    "CS.PA",     # AXA
    "DG.PA",     # Vinci
    "EL.PA",     # EssilorLuxottica
    "SAF.PA",    # Safran
    "KER.PA",    # Kering
    "RI.PA",     # Pernod Ricard
    "BN.PA",     # Danone
    "SGO.PA",    # Saint-Gobain
    "GLE.PA",    # Societe Generale (Often on the edge of inclusion)

    # --- GERMANY (approx. 15) ---
    "SAP.DE",    # SAP SE
    "SIE.DE",    # Siemens AG
    "ALV.DE",    # Allianz SE
    "DTE.DE",    # Deutsche Telekom
    "MBG.DE",    # Mercedes-Benz Group
    "BMW.DE",    # BMW
    "VOW3.DE",   # Volkswagen (Preferred)
    "MUV2.DE",   # Munich Re
    "BAS.DE",    # BASF
    "BAYN.DE",   # Bayer
    "ADS.DE",    # Adidas
    "DHL.DE",    # DHL Group (formerly Deutsche Post DPW.DE)
    "DB1.DE",    # Deutsche Boerse
    "IFX.DE",    # Infineon Technologies
    "ENR.DE",    # Siemens Energy (Recent addition candidate)
    # "DBK.DE",  # Deutsche Bank (Check if re-added in Sep 2025 rebalance)
    
    # --- NETHERLANDS (approx. 6) ---
    "ASML.AS",   # ASML Holding
    "ADYEN.AS",  # Adyen
    "INGA.AS",   # ING Groep
    "AD.AS",     # Ahold Delhaize
    "PHIA.AS",   # Philips
    "PRX.AS",    # Prosus

    # --- SPAIN (approx. 4) ---
    "ITX.MC",    # Inditex
    "IBE.MC",    # Iberdrola
    "SAN.MC",    # Banco Santander (Note: Suffix .MC)
    "BBVA.MC",   # BBVA

    # --- ITALY (approx. 4) ---
    "ISP.MI",    # Intesa Sanpaolo
    "ENEL.MI",   # Enel
    "ENI.MI",    # Eni
    "STLAM.MI",  # Stellantis (Listed in Milan and Paris, .MI often used for index)
    "RACE.MI",   # Ferrari (Recent addition)
    # "UCG.MI",  # UniCredit (Often rotates in/out)

    # --- BELGIUM (1) ---
    "ABI.BR",    # Anheuser-Busch InBev

    # --- FINLAND (2) ---
    "KNEBV.HE",  # Kone
    "NOKIA.HE"   # Nokia
]

# ==========================================
# 2. HELPER FUNCTIONS
# ==========================================
def get_json(endpoint):
    """Helper to handle FMP API requests"""
    url = f"{BASE_URL}/{endpoint}"
    if '?' in endpoint:
        url += f"&apikey={API_KEY}"
    else:
        url += f"?apikey={API_KEY}"
    
    try:
        response = requests.get(url)
        response.raise_for_status()
        return response.json()
    except Exception as e:
        print(f"Error fetching {endpoint}: {e}")
        return []

def get_aggregated_pe(tickers, lookback_years=5):
    """
    Calculates the aggregate P/E of a list of tickers over time.
    (Sum of Market Caps) / (Sum of Net Income)
    """
    print(f"--- Building Index P/E for {len(tickers)} companies ---")
    
    df_earnings = pd.DataFrame()
    df_mcap = pd.DataFrame()
    
    # Limit data to keep script snappy
    limit_q = lookback_years * 4 + 4 
    limit_d = lookback_years * 365
    
    for ticker in tickers:
        # A. Fetch Earnings (Quarterly)
        inc_stmt = get_json(f"income-statement/{ticker}?period=quarter&limit={limit_q}")
        if inc_stmt:
            df = pd.DataFrame(inc_stmt)
            df['date'] = pd.to_datetime(df['date'])
            df.set_index('date', inplace=True)
            df = df.sort_index()
            # Calculate TTM Net Income (Sum last 4 quarters)
            ttm = df['netIncome'].rolling(4).sum()
            # Resample to daily (forward fill)
            df_earnings[ticker] = ttm.resample('D').ffill()

        # B. Fetch Market Cap (Daily)
        mcap = get_json(f"historical-market-capitalization/{ticker}?limit={limit_d}")
        if mcap:
            df = pd.DataFrame(mcap)
            df['date'] = pd.to_datetime(df['date'])
            df.set_index('date', inplace=True)
            df_mcap[ticker] = df['marketCap']

    # C. Aggregate
    # Align dates
    common_idx = df_earnings.index.intersection(df_mcap.index)
    df_earnings = df_earnings.loc[common_idx]
    df_mcap = df_mcap.loc[common_idx]
    
    # Sum components
    total_earnings = df_earnings.sum(axis=1)
    total_mcap = df_mcap.sum(axis=1)
    
    # Calculate P/E
    # Filter out potential data glitches (negative or zero sums)
    idx_pe = total_mcap / total_earnings
    idx_pe = idx_pe[idx_pe > 0] 
    
    return idx_pe


# ==========================================
# 3. MAIN EXECUTION
# ==========================================

# A. Compute STOXX 50 P/E
stoxx_pe = get_aggregated_pe(stoxx_tickers)

# B. (Optional) Get Benchmark. 
# Since S&P 500 P/E is hard to get via free/basic calls, 
# let's generate a synthetic line for comparison based on market averages if API fails,
# or plot just the STOXX history.
# For this script, we will focus on the STOXX Mean Reversion.

# C. Calculate Statistics
current_pe = stoxx_pe.iloc[-1]
avg_pe_3y = stoxx_pe.rolling(window=365*3).mean().iloc[-1]
min_pe = stoxx_pe.min()
max_pe = stoxx_pe.max()

print(f"\nRESULTS:")
print(f"Current STOXX 50 P/E: {current_pe:.2f}")
print(f"3-Year Average P/E:   {avg_pe_3y:.2f}")
print(f"Min P/E (Period):     {min_pe:.2f}")
print(f"Max P/E (Period):     {max_pe:.2f}")

if current_pe < avg_pe_3y:
    print("CONCLUSION: UNDERVALUED relative to recent history.")
else:
    print("CONCLUSION: OVERVALUED relative to recent history.")

# ==========================================
# 4. PLOTTING
# ==========================================
sns.set_style("darkgrid")
plt.figure(figsize=(12, 6))

# Plot Historical P/E
plt.plot(stoxx_pe.index, stoxx_pe, label='STOXX 50 P/E (Calculated)', color='#003399', linewidth=2)

# Plot Average Line
plt.axhline(y=stoxx_pe.mean(), color='r', linestyle='--', label=f'Historical Mean ({stoxx_pe.mean():.1f}x)')

# Add formatting
plt.title('STOXX 50 Historical Valuation (P/E Ratio)', fontsize=16)
plt.ylabel('Price / Earnings Ratio')
plt.xlabel('Year')
plt.legend()

# Highlight Current Valuation
plt.annotate(f'Current: {current_pe:.1f}x', 
             xy=(stoxx_pe.index[-1], current_pe), 
             xytext=(stoxx_pe.index[-1] - pd.Timedelta(weeks=20), current_pe + 2),
             arrowprops=dict(facecolor='black', shrink=0.05))

plt.tight_layout()
plt.show()

And the fact is that it works.

This is impressive to say the least.