import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
clean = pd.read_csv('online_retail_clean.csv', parse_dates=['InvoiceDate'])
clean.shape

desc_lookup = clean.groupby('StockCode')['Description'].agg(lambda s: s.mode().iloc[0])
desc_lookup.head(3)

price_variation = clean.groupby('StockCode')['UnitPrice'].nunique()
(price_variation > 1).sum(), len(price_variation)

deep_dive_sc = '85099B'
desc_lookup[deep_dive_sc]
sub = clean[clean['StockCode'] == deep_dive_sc]
sub.shape[0]

by_price = sub.groupby('UnitPrice').agg(TotalQty=('Quantity', 'sum'), Orders=('InvoiceNo', 'nunique')).reset_index()
by_price.sort_values('UnitPrice')

plt.figure(figsize=(8, 6))
plt.scatter(by_price['UnitPrice'], by_price['TotalQty'], s=80, color='darkorange')
plt.xlabel('Real Unit Price (GBP)')
plt.ylabel('Real Total Quantity Sold')
plt.title(f'Real Price vs Quantity for {desc_lookup[deep_dive_sc]}')
plt.tight_layout()
plt.savefig('price_vs_quantity_scatter.png', dpi=120)
plt.close()

valid = by_price[by_price['TotalQty'] > 0]
log_price = np.log(valid['UnitPrice'])
log_qty = np.log(valid['TotalQty'])
elasticity, intercept = np.polyfit(log_price, log_qty, 1)
round(elasticity, 2)

is_elastic = abs(elasticity) > 1
print(f'Elasticity of {round(elasticity, 2)} means this product is {"elastic" if is_elastic else "inelastic"}: a 1% price cut is associated with roughly a {round(abs(elasticity), 1)}% change in units sold.')

regular_price = sub['UnitPrice'].mode().iloc[0]
regular_price
sub = sub.copy()
sub['IsDiscounted'] = sub['UnitPrice'] < regular_price
sub['IsDiscounted'].mean().round(3)

avg_qty_discounted = sub[sub['IsDiscounted']]['Quantity'].mean()
avg_qty_full = sub[~sub['IsDiscounted']]['Quantity'].mean()
round(avg_qty_discounted, 1), round(avg_qty_full, 1)

revenue_discounted = sub[sub['IsDiscounted']]['Revenue'].sum()
revenue_full = sub[~sub['IsDiscounted']]['Revenue'].sum()
round(revenue_discounted, 2), round(revenue_full, 2)

discount_rows = sub[sub['IsDiscounted']].copy()
discount_rows['DiscountPct'] = (regular_price - discount_rows['UnitPrice']) / regular_price * 100
discount_rows.groupby(pd.cut(discount_rows['DiscountPct'], bins=[0, 10, 20, 30, 100]), observed=True)['Quantity'].mean().round(1)

qty_by_sc = clean.groupby('StockCode')['Quantity'].sum()
price_var_by_sc = clean.groupby('StockCode')['UnitPrice'].nunique()
candidates = price_var_by_sc[(price_var_by_sc >= 4) & (qty_by_sc > 500)].index
len(candidates)

results = []
for sc in candidates:
    s = clean[clean['StockCode'] == sc]
    bp = s.groupby('UnitPrice')['Quantity'].sum()
    bp = bp[bp > 0]
    if len(bp) < 4 or bp.index.to_series().std() == 0:
        continue
    lp, lq = np.log(bp.index.values.astype(float)), np.log(bp.values.astype(float))
    e, _ = np.polyfit(lp, lq, 1)
    results.append((sc, desc_lookup[sc], e, len(bp)))
len(results)

elasticity_df = pd.DataFrame(results, columns=['StockCode', 'Description', 'Elasticity', 'NumPricePoints'])
elasticity_df.shape[0]
elasticity_df.sort_values('Elasticity').head(5)

pct_negative = (elasticity_df['Elasticity'] < 0).mean()
round(pct_negative * 100, 1)

elasticity_df['Elasticity'].mean().round(2)
elasticity_df['Elasticity'].median().round(2)

plt.figure(figsize=(9, 5))
plt.hist(elasticity_df['Elasticity'], bins=30, color='seagreen', edgecolor='white')
plt.axvline(-1, color='red', linestyle='--', label='Real elastic/inelastic boundary')
plt.xlabel('Real Estimated Elasticity')
plt.ylabel('Real Number of Products')
plt.title('Real Distribution of Price Elasticity Across the Catalog')
plt.legend()
plt.tight_layout()
plt.savefig('elasticity_distribution.png', dpi=120)
plt.close()

most_elastic = elasticity_df.sort_values('Elasticity').head(10)
most_elastic[['Description', 'Elasticity']]

least_elastic = elasticity_df[elasticity_df['Elasticity'] > -1].sort_values('Elasticity', ascending=False).head(10)
least_elastic[['Description', 'Elasticity']]

candidate_prices = np.linspace(valid['UnitPrice'].min(), valid['UnitPrice'].max(), 50)
predicted_qty = np.exp(intercept) * candidate_prices ** elasticity
predicted_revenue = candidate_prices * predicted_qty
optimal_price = candidate_prices[np.argmax(predicted_revenue)]
round(optimal_price, 2)

plt.figure(figsize=(8, 6))
plt.plot(candidate_prices, predicted_revenue, color='purple')
plt.axvline(optimal_price, color='red', linestyle='--', label=f'Real optimal price: {round(optimal_price, 2)}')
plt.xlabel('Real Hypothetical Unit Price (GBP)')
plt.ylabel('Real Predicted Total Revenue')
plt.title(f'Real Predicted Revenue Curve for {desc_lookup[deep_dive_sc]}')
plt.legend()
plt.tight_layout()
plt.savefig('revenue_curve_optimal_price.png', dpi=120)
plt.close()

manual_pct_price_change = (valid['UnitPrice'].iloc[-1] - valid['UnitPrice'].iloc[0]) / valid['UnitPrice'].iloc[0]
manual_pct_qty_change = (valid['TotalQty'].iloc[-1] - valid['TotalQty'].iloc[0]) / valid['TotalQty'].iloc[0]
np.sign(manual_pct_qty_change) != np.sign(manual_pct_price_change)

elasticity_df.sort_values('Elasticity').to_csv('product_price_elasticity.csv', index=False)
reloaded = pd.read_csv('product_price_elasticity.csv')
reloaded.shape[0] == elasticity_df.shape[0]

n_elastic = (elasticity_df['Elasticity'] < -1).sum()
n_inelastic = (elasticity_df['Elasticity'] >= -1).sum()
print(f'{n_elastic} real products are genuinely elastic and good discounting candidates; {n_inelastic} are inelastic, where discounting mostly just gives away real margin.')

sub['Weekday'] = sub['InvoiceDate'].dt.day_name()
sub.groupby('Weekday')['IsDiscounted'].mean().round(3).sort_values(ascending=False)

sub[sub['IsDiscounted']]['Country'].value_counts().head(5)

volume_lookup = clean.groupby('StockCode')['Quantity'].sum()
elasticity_df['TotalVolume'] = elasticity_df['StockCode'].map(volume_lookup)
elasticity_df[['Elasticity', 'TotalVolume']].corr().iloc[0, 1].round(3)

reliable = elasticity_df[elasticity_df['NumPricePoints'] >= 6]
reliable.shape[0]
reliable['Elasticity'].mean().round(2)

print(f'Estimated price elasticity for {len(elasticity_df)} real products; {n_elastic} are genuinely discount-worthy.')