# Calculate the mean
df_salary_filtered = df[df['Salary']>10000]

# Group the data by 'Department' and calculate the mean and standard deviation for 'Salary'
department_stats_salary_filtered = df_salary_filtered.groupby('Department').agg({'Salary': ['mean', 'std']})

# Calculate the Coefficient of Variation (COV)
department_stats_salary_filtered['Coefficient of Variation'] = department_stats_salary_filtered['Salary']['std'] / department_stats_salary_filtered['Salary']['mean'] * 100

# Rename the columns to 'Mean' and 'Standard Deviation'
department_stats_salary_filtered.columns = ['Mean', 'Standard Deviation', 'Coefficient of Variation']

# Round the calculated statistics to two decimal places
department_stats_salary_filtered = department_stats_salary_filtered.round(2)

department_stats_salary_filtered.sort_values(by=['Coefficient of Variation'], ascending=False, inplace=False)

department_stats_salary_filtered
print(department_stats_salary_filtered)
# department_stats_salary_filtered.to_csv('comparison.csv', index=True) ## Validated Calulation in SQL query is correct, commented out code block to not save .csv file

# Filter out rows where 'Salary' is < 10000
df_salary_filtered = df[df['Salary'] > 10000]

## Slicing original df to obtain only the Salary column
salary_column_distribution = df_salary_filtered['Salary']
salary_column_distribution

# Calculate the mean and median
mean_salary = round(salary_column_distribution.mean(),3)
median_salary = round(salary_column_distribution.median(),3)

print("Mean Salary:", mean_salary)
print("Median Salary:", median_salary)

# Filter out rows where 'Salary' is < 10000
df_salary_filtered = df[df['Salary'] > 10000]

## Slicing original df to obtain only the Salary column
salary_column_distribution = df_salary_filtered['Salary']

# Calculate sample size
sample_size = df_salary_filtered['Salary'].shape[0]

# Calculate skewness
skewness = df_salary_filtered['Salary'].skew()

# Calculate kurtosis
kurtosis = df_salary_filtered['Salary'].kurtosis()

# # Shapiro-Wilk test for normality, ignore
# shapiro_statistic, shapiro_p_value = stats.shapiro(df_salary_filtered['Salary'])

print("Sample Size:", sample_size)
print("Skewness:", skewness)
print("Kurtosis:", kurtosis)
# print("Shapiro-Wilk Test Statistic:", shapiro_statistic)
# print("Shapiro-Wilk Test p-value:", shapiro_p_value)

# Shapiro-Wilk test for normality
# shapiro_test_statistic, shapiro_p_value = stats.shapiro(salary_column_distribution)
# print(f"Shapiro-Wilk test statistic: {shapiro_test_statistic}, p-value: {shapiro_p_value}")

# Visualization: Histogram
plt.figure(figsize=(8, 5))
sns.histplot(salary_column_distribution, kde=True,
#              hue = 'Departments',
             bins='auto',
#              bins=80
            )
plt.title("Histogram of Annual Salaried Employees")
plt.xlabel("Salary")
plt.ylabel("Frequency")
plt.show()

# Visualization: Q-Q plot
plt.figure(figsize=(8, 5))
stats.probplot(salary_column_distribution, dist="norm", plot=plt)
plt.title("Quantile-Quantile Plot of Annual Salaried Employees")
plt.xlabel("Theoretical Quantiles")
plt.ylabel("Sample Quantiles")
plt.show()

### Old Comments ###
# print("The above Histogram & Quantile-Quantile [Q-Q] plots serves to illustrate the distribution of the data.")
# print("\nExplanation of the tests:")
# print("\nShapiro-Wilk Test:")
# print("\nThe Shapiro-Wilk Test returns a test statistic and a p-value. If the p-value is greater than 0.05 (commonly chosen significance level), it suggests that the data follows a normal distribution.")
# print("\nQuantile-Quantile Plot [Q-Q plot]:")
# print("\nIn the Quantile-Quantile [Q-Q] plot, if the points closely follow the diagonal line, it indicates that the data is close to normal. If the points deviate substantially, it suggests non-normality.")

# print("\nBecause the Z-Test assumes normal distribution of a dataset, and as we can see from both the Q-Q plot & histogram plot, that salary does not follow a normal distribution across the entire dataset. I've decided to opt for ")

# print("\nDue to the following reasons below, I'll be opting for outlier detection via IQR/boxplot method instead for salary disparity within departments:")
# print("\n 1)  Z-Test assumes normal distribution of a dataset, and as we can see from both the Q-Q plot & histogram plot, that salary does not follow a normal distribution across the entire dataset. Most salary histograms. Most salary plots, including the one obtained from this dataset, is usually right-skewed. This is because most people earn in the low/medium range of salaries, with a few exceptions.")
# print("\n 2) Further investigation revealed that 1155 rows have a salary of <$1000, which is the reason for the left peak in the Histogram plot. This is a large reason for the dataset to not follow a normal distribution [?]")
# print("\n 3) p-values may not be accurate for N > 5000. Our dataset contains 6946 rows in salary")

# Generation of an example of a normally distributed dataset
np.random.seed(3)
mean = 0
std_dev = 1
sample_size = 1000
dummy_data = np.random.normal(loc=mean, scale=std_dev, size=sample_size)

# Shapiro-Wilk test for normality
statistic, p_value = stats.shapiro(dummy_data)
alpha = 0.05

# Create side-by-side subplots for histogram and Q-Q plot
plt.figure(figsize=(15, 6))

# Plot the histogram (normalized)
plt.subplot(1, 2, 1)
plt.title("Histogram of a Typical Normally Distributed Dataset")
sns.histplot(dummy_data, kde=True, color='skyblue', bins=30)
plt.xlabel("Value")
plt.ylabel("Frequency")

# Visualization: Q-Q plot
plt.subplot(1, 2, 2)
stats.probplot(dummy_data, dist="norm", plot=plt)
plt.title("Quantile-Quantile Plot of a Typical Normally Distributed Dataset")
plt.xlabel("Theoretical Quantiles")
plt.ylabel("Sample Quantiles")

# Display Shapiro-Wilk test result as main title
# if p_value > alpha:
#     plt.suptitle(f"Shapiro-Wilk Test: p-value = {p_value:.4f}\nData appears to be normally distributed", y=1.02)
# else:
#     plt.suptitle(f"Shapiro-Wilk Test: p-value = {p_value:.4f}\nData does not appear to be normally distributed", y=1.02)

plt.tight_layout()
plt.show()

import matplotlib.pyplot as plt
import seaborn as sns

df = pd.read_csv('data/Employee_Salaries.csv')
# Filter out rows where 'Salary' is <= 10000
df_salary_filtered = df[df['Salary'] > 10000]
df_salary_filtered['Department'].unique() ## checking departments

# # # Define the list of departments to analyze, in this case I'm inspecting all histograms for annual salaried workers
department_list_to_analyse = ['AGR', 'AUD', 'CAD', 'CCC', 'CIR', 'CIT', 'CLK', 'CMD', 'COM',
       'COR', 'CUL', 'CVB', 'CWA', 'ECC', 'ECO', 'EMS', 'FIN', 'FIR',
       'GRD', 'HNP', 'HRD', 'HSD', 'JUV', 'LIB', 'MSB', 'MUS', 'OEM',
       'PAR', 'PHD', 'PLN', 'POL', 'PUD', 'PWD', 'REA', 'RMO', 'SHF',
       'STR', 'TRE']

# Filter the DataFrame to include only the departments in department_list_to_analyse
df_salary_filtered_dept_slice = df_salary_filtered[df_salary_filtered['Department'].isin(department_list_to_analyse)]

# Define a color palette with unique colors for each department
color_palette = sns.color_palette("Set1", n_colors=len(department_list_to_analyse))

# Loop through each department in department_list_to_analyse and plot a histogram
for i, department in enumerate(department_list_to_analyse):
    plt.figure(figsize=(8, 5))
    sns.histplot(data=df_salary_filtered_dept_slice[df_salary_filtered_dept_slice['Department'] == department], 
                 x='Salary', kde=True, bins='auto', color=color_palette[i], legend=True)
    plt.title(f"Histogram of Salary for Department: {department}")
    plt.xlabel("Salary")
    plt.ylabel("Frequency")
#     plt.savefig(f"data/exported/AnnualSalary_Histogram_{department}.png") # Saves each department as a separate .png file, commented out for presentation
    plt.show()

df_Stats = pd.read_csv('data/EmployeeSalaries_stddev_Avg_CV_Outlier.csv')
df_Stats

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

# Display basic information about the DataFrame
print(df_Stats.info())

# Check for missing values
print(df_Stats.isnull().sum())

# Summarize numerical columns
print(df_Stats.describe())

# Filter out rows where 'Salary' is > 10000
df_salary_filtered = df[df['Salary'] > 10000]

department_list_to_analyse = ['CMD','AGR','CCC','PWD','CAD']
df_salary_filtered_dept_slice = df_salary_filtered[df_salary_filtered['Department'].isin(department_list_to_analyse)]
df_salary_filtered_dept_slice['Department'].unique() ## checking if slice is correct

## Visualization: Histogram for each department
plt.figure(figsize=(12, 6))
sns.histplot(data=df_salary_filtered_dept_slice, x='Salary', kde=True, 
             hue='Department',
             bins='auto')
plt.title("Histogram of Annual Salaried Employees by Department")
plt.xlabel("Salary")
plt.ylabel("Frequency")
plt.show()

## plotting salary distribution within departmental level, we also see that they are mostly right-skewed as well for the departments chosen
import matplotlib.pyplot as plt
import seaborn as sns

# Filter the DataFrame to include only salaries greater than or equal to 10000
df_salary_filtered = df[df['Salary'] > 10000]

# Define the list of departments to analyze
department_list_to_analyse = ['CMD', 'AGR', 'CCC', 'PWD', 'CAD']

# Filter the DataFrame to include only the departments in department_list_to_analyse
df_salary_filtered_dept_slice = df_salary_filtered[df_salary_filtered['Department'].isin(department_list_to_analyse)]

# Define a color palette with unique colors for each department
color_palette = sns.color_palette("Set1", n_colors=len(department_list_to_analyse))

# Loop through each department in department_list_to_analyse and plot a histogram
for i, department in enumerate(department_list_to_analyse):
    plt.figure(figsize=(8, 5))
    sns.histplot(data=df_salary_filtered_dept_slice[df_salary_filtered_dept_slice['Department'] == department], 
                 x='Salary', kde=True, bins='auto', color=color_palette[i], legend=True)
    plt.title(f"Histogram of Salary for Department: {department}")
    plt.xlabel("Salary")
    plt.ylabel("Frequency")
    plt.show()

# Box plot of overall dataset
sns.boxplot(data=df)

# Show the plot
plt.show()

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

## iqr_scale_value set to 1.5. Assuming a normally distributed dataset having a Gaussian Distribution, having an iqr_scale_value of 1.5 would theoretically capture outliers that fall outside the +/- 3 Std. Dev from the mean. This is commonly used in outlier detection with a typical Gaussian Distribution, which would typically represent the 0.15% of outlying data outside of a +/- 3 Std. Dev of a normally distributed graph.
## Further tuning of the optimal iqr_scale_value to accurately decpict outlier detection can be explores, given this dataset is not normally distributed
IQR_SCALE_VALUE = 1.5 ## parameter that controls the sensitivity of outlier detection, a higher value will result in fewer outliers being counted from both upper_bound & lower_bound parameters

# Calculate the number of outliers for each department
def count_outliers(s):
    q1 = s.quantile(0.25)
    q3 = s.quantile(0.75)
    iqr = q3 - q1
    lower_bound = q1 - IQR_SCALE_VALUE * iqr ## Define lower threshold limit as lower_bound, graphically this will affect the lower whisker length of the boxplot
    upper_bound = q3 + IQR_SCALE_VALUE * iqr ## Define upper threshold limit as upper_bound, graphically this will affect the upper whisker length of the boxplot
    return s[(s < lower_bound) | (s > upper_bound)].count() ## any point that falls outside the lower_bound and upper_bound will be treated as an 'outlier' point and be counted

outliers_count_by_department = df.groupby('Department')['Salary'].apply(count_outliers)
print(f"Count of outliers by departments in descending order via IQR-method:\n {outliers_count_by_department.sort_values(ascending=False)}") ## printout of departments by DESC order based on the count of outliers beyond the upper threshold and 

# Sort the departments based on the count of outliers in descending order
sorted_departments = outliers_count_by_department.sort_values(ascending=False).index

#### Filter out the outliers based on iqr_scale_value ####
Q1 = df['Salary'].quantile(0.25)
Q3 = df['Salary'].quantile(0.75)
IQR = Q3 - Q1
df_filtered = df[~((df['Salary'] < (Q1 - IQR_SCALE_VALUE * IQR)) | (df['Salary'] > (Q3 + IQR_SCALE_VALUE * IQR)))]

# Create the box plot segmented by 'Department' in the sorted order using the filtered data
plt.figure(figsize=(20, 8))
sns.boxplot(x='Department', y='Salary', data=df_filtered, order=sorted_departments)
plt.title("IQR Plot: Salary Distribution by Department (Descending Order of Outliers)")
plt.xlabel("Department")
plt.ylabel("Salary")

# ### ORIGINAL BOX PLOT without taking into account of iqr_scale_value ####
# plt.figure(figsize=(20, 8))
# sns.boxplot(x='Department', y='Salary', data=df, order=sorted_departments)
# ## Add a title and labels
# plt.title("IQR Plot: Salary Distribution by Department (Descending Order of Outliers)")
# plt.xlabel("Department")
# plt.ylabel("Salary")

### Rotate the x-axis labels for better visibility
### plt.xticks(rotation=45, ha='right')

#### Show the plot
plt.show()