Applying FAMD in Python: Preparing Mixed Data for Clustering

From data exploration to dimensionality reduction using the Prince library – building the foundation for K-Means clustering.

Photo by Bioscience Image Library by Fayette Reynolds on Unsplash

Introduction

In the first part of this series (From PCA to FAMD: Dimensionality Reduction for Mixed Data | by Georgios Kokkinopoulos | Apr, 2026 | DataDrivenInvestor), we explored the intuition behind PCA and why Factor Analysis of Mixed Data (FAMD) is better suited for datasets containing both numerical and categorical variables.

In this article, we move from theory to practice. Using a real-world customer dataset, we will prepare mixed data for clustering by applying FAMD in Python.

We begin with data exploration, handle missing values, and then apply FAMD using the Prince library to transform the dataset into a lower-dimensional space.

Finally, we use a scree plot to determine the optimal number of components — creating a clean and meaningful representation of the data that will be used for clustering in the next part of the project.

By the end of this article, you will have a transformed dataset ready for clustering.

What dataset will be used

The dataset we will be using was initially created for a competition in Customer Segmentation in Analytics Vidhya site Customer Segmentation Hackathon and was retrieved by Kaggle which is a top destination for Data Scientists and Machine Learning practitioners looking to participate in a competition or to find a robust dataset to build a model on. The dataset was stored there by user Vetrivel-PS | Kaggle. Link to the dataset here: Customer Segmentation. It contains demographic (Age, Gender etc.) and behavioural (Spending Score) data of 8068 customers of an automotive company, so it effectively reflects the nature of mixed data that marketing people work with in real-life cases.

Data exploration

Read dataset — Identify missing values

# import libraries
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

# Read the customers dataset
df=pd.read_csv('customers.csv')
df.head()

print('Number of rows and columns: '+str((df.shape)))
print('')
print('Column types:')
print('')
print(df.info())

Missing values identified in some columns so this will be examined in the features’ visualisation.

Exploring each feature

Let’s define three functions that will be used for visualisation:

def aggr_table(df,x,sort=False):
    """Creates count,cumulative count and the equivalent % columns for a specific field in a df

    Args:
    df: The dataframe
    x: string. The column in the dataframe we group by

    Returns:
    A table with count, cumulative count, count% and cumulative count% for the given column.
    """

    aggr=df.groupby(x).size().reset_index(name='Count')    # Create the groupby and the first metric
    y_per_cent='Count %'    # Will be used as the name of the fourth metric
    aggr[y_per_cent]=round(100*aggr['Count']/sum(aggr['Count']),1)    # Count %

    if not sort and pd.api.types.is_numeric_dtype(df[x]):
        cum_name='Cumulative Count'   # Name of the second metric created
        cum_name_per_cent=cum_name+' %'    # Name of the third metric created
        aggr[cum_name]=aggr['Count'].cumsum()    # Cumulative count created
        aggr[cum_name_per_cent]=round(100*aggr[cum_name]/sum(aggr['Count']),1)    # Cumulative count %


    return aggr

def metric_convert(metric):
    """This function will be used for the function below. It just saves the user from submitting long argument names

    Args:
    metric: string.The argument that will be passed in the function aggr_graph below. 

    Returns:
    The full name of the metric
    """

    if metric=='Cum Count':
        metric='Cumulative Count'

    elif metric=='Cum Count %':
        metric='Cumulative Count %'

    elif metric=='Count' or metric=='Count %':
        metric=metric

    else: 
        print('Please give a valid metric name')
        sys.exit()

    return metric     

def aggr_graph(df, x, metric, kind, show_val=True, h_line=None, bins=5, return_table=False, sort=False, ascend=True,missing_mode='detailed'):

    """Creates a graph for a specific feature of a dataframe following a specific metric

    Args:
    df: The dataframe
    x: The feature being visualised (string)
    metric: The metric used in the visualisation. Accepted values are the below ONLY:
            'Count','Count %','Cum Count' for Cumulative Count and 'Cum Count %' for Cumulative Count %
    kind: The kind of plot. Please give one of these ONLY: 'pie','bar','barh' or 'hist'.
    bins: Number of categories for the histogram
    show_val: If True, the values of the feature will appear in the graph. 
    h_line: Optional. Creates horizontal line
    return_table: If True it will also provide the table of the values in the graph
    sort: For bar chart ONLY: If True it will sort the bars by value rather than alphabetically.
    ascend: Use it only if sort argument is True. It's the order by which bars will be sorted. Use True for ascending
    missing_mode: Used in pie chart only. Detailed: Adds missing as a new category. Binary: Presents missing vs non-missing 

    Returns:
    The graph. Pie chart, bar chart or histogram.
    """


    aggr=aggr_table(df,x,sort)    # The aggregate by this feature created. 
    if sort:
        aggr=aggr.sort_values(by='Count',ascending=ascend)    # Aggr table sorted if needed


    if kind == 'pie':
        total = len(df)
        missing_count = df[x].isna().sum()
        non_missing_count = total - missing_count

        if missing_mode == 'binary':
            # Binary pie: Missing vs Non-missing
            pie_df = pd.DataFrame({
                x: ['Non-missing', 'Missing'],
                'Count': [non_missing_count, missing_count],
                'Count %': [
                    non_missing_count / total * 100,
                    missing_count / total * 100
                ]
            })

        else:
            # Detailed pie: categories + Missing
            pie_df = aggr_table(df, x)

            if missing_count > 0:
                missing_row = {
                    x: 'Missing',
                    'Count': missing_count,
                    'Count %': missing_count / total * 100
                }
                pie_df = pd.concat(
                    [pie_df, pd.DataFrame([missing_row])],
                    ignore_index=True
                )

        pie_df.plot(
            kind='pie',
            y='Count %',
            labels=pie_df[x],
            autopct='%1.1f%%',
            legend=False
        )

        title = f'Distribution of {x}'

        plt.ylabel('')  # Removes the default y-axis label for clarity

      # Move legend outside the pie chart
        plt.legend(title=x, bbox_to_anchor=(1.05, 1), loc='upper left')
        plt.tight_layout()

    elif kind=='bar' or kind=='barh':
        metric=metric_convert(metric)    # Function above used to convert the metric to its full name

        plt.figure(figsize=(12, 12))
        aggr.plot(x=x,y=metric,kind=kind)    # Graph created

        plt.xticks(rotation=45)

        # Labels created below

        if kind == 'bar':
            plt.xlabel(x)
            plt.ylabel(metric)
        else:
            plt.xlabel(metric)
            plt.ylabel(x)

        title='Count of Customers by '+x  

    elif kind == 'hist':


        fig, ax = plt.subplots(figsize=(8, 6))

        data = df[x].dropna()

        min_val = data.min()
        max_val = data.max()

        bin_edges = np.linspace(min_val - 0.5,max_val + 0.5,bins + 1)


        counts, bin_edges, patches = ax.hist(data,bins=bin_edges,edgecolor='black')


        ax.set_xlabel(x)
        ax.set_ylabel('Frequency')
        title = f'Count of Customers by {x}'

        # --- BIN LABELS ---
        bin_centers = (bin_edges[:-1] + bin_edges[1:]) / 2

        bin_labels = [f"{int(np.ceil(bin_edges[i]))}–{int(np.floor(bin_edges[i+1]))}"
            for i in range(len(bin_edges) - 1)]


        ax.set_xticks(bin_centers)
        ax.set_xticklabels(bin_labels, rotation=45, ha='right')

        # --- BIN FREQUENCIES ---
        if show_val:
            for count, patch in zip(counts, patches):
                if count > 0:
                    x_pos = patch.get_x() + patch.get_width() / 2
                    ax.text(
                        x_pos,
                        count,
                        int(count),
                        ha='center',
                        va='bottom',
                        fontsize=9
                    )



    # If statements below create horizontal line and show the values if requested 
    if h_line is not None:
        plt.axhline(y=h_line, color='r', linestyle='--')
    if show_val and kind == 'barh':
        for i, value in enumerate(aggr[metric]):
           plt.text(value * 0.95, i, str(value), va='center', ha='right', color='black')

    elif show_val and kind!='pie' and kind!='hist':
        for i, value in enumerate(aggr[metric]):
           plt.text(i, value , str(value), ha='center', va='bottom')

    plt.title(title)    # Title added

    plt.show()    # Display the chart

    # Table with values appears below the graph if requested
    if return_table:
        return aggr

Let’s start visualising the features

# Distribution of Age
aggr_graph(df,'Age','Count','hist',bins=5)

Positively skewed distribution — Median is at around 40 years of age

# Distribution of Gender
aggr_graph(df,'Gender','Count %','pie')

More Male than Female customers without the former being dominant though.

# Distribution of Ever Married
aggr_graph(df,'Ever_Married','Count %','pie')

Two distinct categories here (not a very dominant one) plus a substantial proportion of missing values. That could be a third category in our final segmentation.

# Distribution of Graduated or not
aggr_graph(df,'Graduated','Count %','pie')

More Graduated than non-Graduated customers, but the latter form a very considerable category too.

# Distribution of Profession
aggr_graph(df,'Profession','Count','bar',sort=True,ascend=False)

Artist is the biggest category here (almost double compared to the second one), but even the least populated categories (Marketing, Homemaker) should be included in the model.

It’s also worth examining what proportion of customers didn’t provide information about their profession. Let’s see the binary pie-chart below:

aggr_graph(df,'Profession','Count %','pie',missing_mode='binary')

Another 1.5% missing values for this feature — not to be ignored

# Distribution of Work Experience
aggr_graph(df,'Work_Experience','Count','hist',bins=5,return_table=True)

The vast majority of Work Experience values are between 0 and 2 years. Let’s also explore the missing values:

# Looking for missing values in Work Experience
aggr_graph(df,'Work_Experience','Count %','pie',missing_mode='binary')

Missing is a very considerable proportion of customers. This will definitely have to be imputed.

# Distribution of Spending Score
aggr_graph(df,'Spending_Score','Count','barh',sort=True,ascend=True)

All 3 categories are considerably densely populated.

# Distribution of Family Size
aggr_graph(df,'Family_Size','Count','bar')

Positively skewed with a median at around 3

aggr_graph(df,'Family_Size','Count %','pie',missing_mode='binary')

Missing is again a remarkable proportion which can’t be dropped from the dataset.

Data pre-processing — Handling missing values

As we saw above, there are 5 variables with missing values: Two numeric ones (Work Experience and Family Size) and three categorical ones (Ever Married, Graduated and Profession). We are going to handle these two types in different ways:

• Numeric features: Median imputation. By filling missing values with the median instead of the mean we are achieving two benefits:

1. We preserve the skewness of the distribution. By assigning missing values to the mean we would shift the distribution towards the direction of the skewness
2. We do not shrink the variance too aggressively. Missing values filled by the mean would mean that all new points fall on the average resulting in no contribution to the variance at all. Remember that variance is needed in clustering as this will help us distinguish between clusters.

• Categorical features: We will treat “missing” as a new category. There are two main reasons behind this:

1. We preserve uncertainty. By assigning a missing value to one of the existing values based on the values of other features (e.g. the majority of people over 70 are Lawyers so we will assign the missing job of any individual over 70 to the value of Lawyer) we artificially increase the homogeneity of a cluster which mainly consists of Law professionals in their 70s. That should be avoided since there is no safe way we can conclude that those who didn’t provide this piece of information have the same profession as other people with similar characteristics.
2. Missingness is often a piece of information about these individuals: Those who prefer not to share such details. So category “missing” should be considered

Below is the code we are using to process these features:

# Create a copy of the initial dataset to work on
df_clust=df.copy()

# Fill missing values in categorical variables with a new category ("Missing") 
# Change type from Object to Category

cat_cols = df_clust.select_dtypes(include='object').columns

for col in cat_cols:
    df_clust[col] = df_clust[col].fillna("Missing").astype("category")

# Fill missing values of numeric variables with the Median 

num_cols = df_clust.select_dtypes(include='number').columns

for col in num_cols:
    if df_clust[col].isnull().sum()>0:

        median = df_clust[col].median()
        df_clust[col] = df_clust[col].fillna(median).astype('float')

# Drop ID as not needed for clustering
df_clust=df_clust.drop('ID',axis=1)

Why Dimension Reduction

Before actually applying FAMD to the original dataset it’s worth mentioning some of the benefits of Dimension Reduction in K-Means and other clustering algorithms, especially when we have mixed data with numerous categories:

• Makes the algorithm computationally efficient — especially useful for large datasets.
• Reduces the noise from the transformed dataset while it captures the structure hidden in it by focusing on the components that account for the most variance.
• Provides the opportunity for 2-D and 3-D visualisation of the dataset, subsequently of the clusters by using the first two or three Principal Components respectively. You will see a 2-D representation in the next paragraph.
• Mitigates the “Curse of Dimensionality” which would otherwise result in sparse data points that prevent the algorithm from identifying compact clusters.

Dimension Reduction implementation

FAMD — Method implementation

Let’s import our libraries first. We will be using the library “prince” for the FAMD algorithm as this is not available in SKlearn. Then, we will initiate FAMD by using all 8 initial variables as the initial components. FAMD handles directly the initial variables (hence n_components=8) and does the one-hot encoding for us. Finally, we will apply the FAMD transformation to the clean dataset we have created for our project.

# Import libraries
import prince

# FAMD: Factor Analysis for Mixed Data 
# We have 8 variables in total so we use 8 components
famd = prince.FAMD(n_components=8, random_state=42)
df_famd = famd.fit_transform(df_clust)

Now we have created the transformed dataset, let’s take a look at the scatterplot of the first two Principal Components. This should provide an initial view of the potential clusters that will be created later.

# Create a new figure with a specific size
plt.figure(figsize=(10, 6))

# Extract the first two FAMD components
comp1 = df_famd.to_numpy()[:, 0]
comp2 = df_famd.to_numpy()[:, 1]

# Create a scatter plot for an initial view of potential clusters
scatter = plt.scatter(comp1, comp2,alpha=0.2)  
plt.title('Transformed Dataset Visualisation')   #Add title

# Label axes with explained variance
plt.xlabel(f'Comp 1 ({explained_inertia[0]*100:.1f}% explained)')
plt.ylabel(f'Comp 2 ({explained_inertia[1]*100:.1f}% explained)')

# Add grid
plt.grid(True)

plt.show()

Image by Author

It appears that we will finally create at least two densely populated clusters at the bottom and one or more scarcely populated clusters at the top.

How to determine the number of principal components

The first step is to calculate the explained inertia of each of the 8 principal components that the FAMD has created. This is the proportion of the total variance in the dataset explained by each of the PCs which by definition of dimension reduction (PCA or FAMD) will be in descending order. The code below addresses it.

# Extract eigenvalues from the fitted FAMD model
# Each eigenvalue represents the variance explained by a component
eigenvalues = famd.eigenvalues_

# Convert eigenvalues to explained inertia (proportion of total variance)
# This shows how much information each component captures
explained_inertia = eigenvalues / eigenvalues.sum()

# Round explained inertia values to 2 decimal places for readability
explained_inertia_2dp = [round(x, 2) for x in explained_inertia]

# Display explained inertia per component
print('Explained Inertia:')
print(explained_inertia_2dp)
print('')

# Compute cumulative explained inertia across components
# This helps determine how many components are needed to capture most variance
explained_inertia_cum_2dp = [round(x, 2) for x in np.cumsum(explained_inertia)]

# Display cumulative explained inertia
print('Cumulative Explained Inertia:')
print(explained_inertia_cum_2dp)

The explained inertia above gives us some insight as to how many of the 8 PCs we should use in our K-Means algorithm. The general rule is : Keep adding components as long as the last one included adds a significant amount of inertia to the model. A bar chart (scree plot) will definitely be a better guide for us than mere numbers. The code below will provide it for us:

# Create component index for the scree plot
components = range(famd.n_components)

# Plot bar chart of explained inertia per component
bars = plt.bar(components, explained_inertia, color="blue")

# Add value labels on top of each bar
for bar, val in zip(bars, explained_inertia):
    plt.text(
        bar.get_x() + bar.get_width() / 2,  # Center of the bar
        bar.get_height(),                   # Height of the bar
        f"{val:.2f}",                       # Format value to 2 decimals
        ha='center',
        va='bottom'
    )

# Axis labels
plt.xlabel('FAMD components')
plt.ylabel('Inertia')

# Set x-axis ticks
plt.xticks(components)

# Display the plot
plt.show()

Scree plot below:

The bar chart gets flat from the 4th component on. That’s an indication that we need to keep 4 components only as any additional ones will probably add more noise than explained variance. Total variance explained by these 4 components is 64% (see Cumulative Explained Inertia above) which is satisfactorily enough as it says that these 4 components will explain roughly two thirds of the total variance.

Conclusion

In this article we walked through the following steps:

• Initial visualisation and imputation of missing values (median for numeric — standalone category for categorical variables).
• Dimension Reduction implementation by the use of the FAMD method.
• How to determine the number of Principal Components that we will keep from FAMD in order to create the transformed dataset.

In the third part of this analysis we will apply K-Means in the transformed dataset to achieve our goal: Interpretable and meaningful for the business clusters.

If you found this article useful please give me a like and don’t hesitate to follow me. You are also more than welcome to connect with me on LinkedIn (Georgios Kokkinopoulos | LinkedIn).

References

[1] William Blaufuks, FAMD: How to generalize PCA to categorical and numerical data | by William Blaufuks | TDS Archive | Medium (2021), TDS Archive

[2] Calvin Hui, The Math Behind Principal Component Analysis (PCA): Variance, Reconstruction Error, Eigenvectors, and Intuition | Data Science Collective (2025), Data Science Collective

[3] Dario Radecic, 5 PCA Visualizations You Must Try On Your Next Data Science Project | by Dario Radečić | TDS Archive | Medium (2024), TDS Archive

[4] Dataset: Customer Segmentation CC0:Public Domain

[5] Github repo: GKokkinopoulos/Customer-Segmentation