#######################
# Logistic Regression #
#######################

# This script reproduces the full Chestnut Ridge logistic regression example in Python:
#   Step 1. Logistic regression model (fit, pseudo R-square, odds ratios, predicted probabilities)
#   Step 2. Breakeven rate for the catalog campaign
#   Step 3. Predicted vs. actual purchases
#   Step 4. Confusion matrix and hit rate
#   Step 5. Profitability and ROMI compared to random and RFM targeting
# Run it from top to bottom. Each step prints its result or shows a plot.
# The "# %%" lines split the script into cells you can run one at a time
# in VS Code, Spyder, or PyCharm.

# %%
## Install Packages (if needed - only installs packages that are missing)
## Each entry is {name used with import: name used with pip install}
import sys, subprocess, importlib.util
packages = {"pandas": "pandas", "numpy": "numpy", "statsmodels": "statsmodels", "matplotlib": "matplotlib"}
new_packages = [pip_name for import_name, pip_name in packages.items() if importlib.util.find_spec(import_name) is None]
if new_packages: subprocess.check_call([sys.executable, "-m", "pip", "install", *new_packages])

## Load Packages and Set Seed
from tkinter import Tk, filedialog
import numpy as np
import pandas as pd
import statsmodels.api as sm
import statsmodels.formula.api as smf
import matplotlib.pyplot as plt

np.random.seed(1)

## Functions to choose a file to open, or name a file to save, with a dialog box
## (the Python version of file.choose() in R - no need to set the working directory)
def file_choose(title="Choose a file"):
	root = Tk()
	root.withdraw()                     # hide the empty Tk window
	root.attributes("-topmost", True)   # bring the dialog to the front
	path = filedialog.askopenfilename(title=title, filetypes=[("CSV files", "*.csv"), ("All files", "*.*")])
	root.destroy()
	if not path:
		raise SystemExit("No file was chosen.")
	return path

def file_save(title="Save the file as", initialfile=""):
	root = Tk()
	root.withdraw()
	root.attributes("-topmost", True)
	path = filedialog.asksaveasfilename(title=title, initialfile=initialfile, defaultextension=".csv",
		filetypes=[("CSV files", "*.csv")])
	root.destroy()
	if not path:
		raise SystemExit("No file name was chosen.")
	return path

## Read in Logistic Regression data
logit = pd.read_csv(file_choose("Choose retail_logit.csv")) ## Choose retail_logit.csv file


# %%
############################################
# Step 1: Logistic Regression Model        #
############################################

## Transform and Create Data
logit["number_of_orders2"] = logit["number_of_orders"] ** 2
logit["lnrevenue"] = np.log(logit["revenue"] + 1)

## Run Logistic Regression using GLM
logit_result = smf.glm(formula="purchase ~ lnrevenue + number_of_orders + number_of_orders2 + recency_days + "
	"loyalty_card + married + income", data=logit, family=sm.families.Binomial()).fit()
print(logit_result.summary())

## Pseudo R-square - McFadden
null_result = smf.glm(formula="purchase ~ 1", data=logit, family=sm.families.Binomial()).fit()
print(1 - logit_result.llf / null_result.llf)

## Odds Ratio
print(np.exp(logit_result.params))

## Predicted Probability
logit["predict"] = logit_result.predict(logit)


# %%
############################################
# Step 2: Breakeven Rate                   #
############################################

## Campaign Economics
avg_purchase = 40     ## Average purchase amount
cogs_rate = 0.48      ## Cost of goods sold as a share of revenue ($19.20 / $40)
ship_cost = 6         ## Average cost to ship the product
mail_cost = 2         ## Average cost of the marketing campaign (catalog) per customer

## Profit per Purchasing Customer
profit_per_buyer = avg_purchase - avg_purchase * cogs_rate - ship_cost - mail_cost
print(profit_per_buyer)

## Breakeven Rate
breakeven = mail_cost / profit_per_buyer
print(breakeven)


# %%
############################################
# Step 3: Predicted vs. Actual             #
############################################

## Target Customers at or above the Breakeven Rate
logit["predict_purchase"] = np.where(logit["predict"] >= breakeven, 1, 0)

## Predicted vs. Actual Purchase for Each Customer
print(logit[["customer_id", "predict", "predict_purchase", "purchase"]].head(20))

## Predicted vs. Actual Purchase, Sorted by Predicted Probability
logit_sorted = logit.sort_values("predict", ascending=False)
print(logit_sorted[["customer_id", "predict", "predict_purchase", "purchase"]].head(20))

## Distribution of Predicted Probability for Purchasers and Non-Purchasers
bins = np.arange(0, 1.025, 0.025)
fig, ax = plt.subplots(figsize=(9, 6))
ax.hist(logit.loc[logit["purchase"] == 0, "predict"], bins=bins, alpha=0.6,
	color="#3366CC", edgecolor="white", label="Did Not Purchase")
ax.hist(logit.loc[logit["purchase"] == 1, "predict"], bins=bins, alpha=0.6,
	color="#E68A1A", edgecolor="white", label="Purchased")
ax.axvline(breakeven, color="black", linestyle="--", linewidth=2)
ax.text(breakeven + 0.01, ax.get_ylim()[1] * 0.95,
	f"Breakeven Rate = {np.floor(breakeven * 10000 + 0.5) / 100:.2f}%", va="top")
ax.set_title("Predicted Probability of Purchase by Actual Purchase", loc="left")
ax.set_xlabel("Predicted Probability of Purchase")
ax.set_ylabel("Number of Customers")
ax.legend(loc="upper right", frameon=False)
ax.grid(color="#EBEBEB")
ax.set_axisbelow(True)
for side in ["top", "right", "bottom", "left"]:
	ax.spines[side].set_visible(False)
plt.tight_layout()
plt.show()


# %%
############################################
# Step 4: Confusion Matrix                 #
############################################

## Confusion Matrix (Hit Rate Table)
confusion = pd.crosstab(logit["predict_purchase"], logit["purchase"],
	rownames=["Predicted"], colnames=["Actual"])
print(confusion)

true_negative = confusion.loc[0, 0]
true_positive = confusion.loc[1, 1]
false_negative = confusion.loc[0, 1]
false_positive = confusion.loc[1, 0]

## Hit Rate
hit_rate = (true_negative + true_positive) / confusion.values.sum()
print(hit_rate)


# %%
############################################
# Step 5: Profitability and ROMI           #
############################################

## Function to Compute Campaign Profitability
def campaign_results(n_customers, n_targeted, n_buyers):
	revenue = n_buyers * avg_purchase
	cogs = revenue * cogs_rate
	shipping = n_buyers * ship_cost
	marketing = n_targeted * mail_cost
	profit = revenue - cogs - shipping - marketing
	return pd.Series({"number_of_customers": n_customers,
		"number_targeted": n_targeted,
		"pct_targeted": n_targeted / n_customers,
		"number_of_buyers": n_buyers,
		"response_rate": n_buyers / n_targeted,
		"revenue": revenue,
		"cogs": cogs,
		"shipping": shipping,
		"targeted_marketing": marketing,
		"profit": profit,
		"romi": profit / marketing})

## Random Targeting (send the catalog to every customer)
random = campaign_results(len(logit), len(logit), logit["purchase"].sum())

## Targeted: Independent Sort RFM (results from Chapter 7)
rfm_independent = campaign_results(10000, 4459, 1314)

## Targeted: Sequential Sort RFM (results from Chapter 7)
rfm_sequential = campaign_results(10000, 4480, 1306)

## Targeted: Logistic Regression
logistic = campaign_results(len(logit), true_positive + false_positive, true_positive)

## Compare Results
results = pd.DataFrame({"random": random, "rfm_independent": rfm_independent,
	"rfm_sequential": rfm_sequential, "logistic": logistic})
print(results.round(4))

## Export Logistic Regression Results (used as the input data for the Chapter 9 example)
logit.to_csv(file_save("Save as logit_pred.csv", "logit_pred.csv"), index=False) ## Name file logit_pred.csv