I have shown you an example of a central composite design in the previous blog post: A step by step example of a Central Composite Design (CCD). In this post I am going to show you how you can create the central composite design we used in the example with Python. We use the pyDOE3 package to do that. You can find the documentation here.

Step 1: Installation

You need pyDOE3 to generate the CCD. So if you haven’t already, you need to install it.

pip install pyDOE3

Step 2: Generate a CCD in coded units

We’ll create a CCD for the four factors from the filtration example in our CCD overview:

• T — Temperature: 20–40 °C
• P — Pressure: 1–3 bar
• CoF — Formaldehyde concentration: 2–6 %
• RPM — Agitation speed: 100–300

Two common CCD variants:

• Face-centered (α = 1.0). Good when extreme star points are impractical.
• Rotatable (α = 2^(k/4)). Gives uniform prediction variance at equal radius.

The snippet below defaults to a face-centered CCD with 5 center points. Switch to rotatable by uncommenting the line indicated. Be aware that for a rotatable design the axial points extend outside the factorial space. So you need to make sure that this is actually possible (see next section). For example in our example, a rotatable design (α = 2.0) would require testing at 10 °C, 0 bar, 0 % formaldehyde, and 0 rpm, which is not feasible in practice.

# Build a CCD in coded units
from pyDOE3 import ccdesign
import pandas as pd
import numpy as np

# 1) Define factors and their real-unit ranges
#    These will be used later to map coded → real values.
factors = {
    "T_C": (20, 40),      # Temperature [°C]
    "P_bar": (1, 3),      # Pressure [bar]
    "CoF_pct": (2, 6),    # Formaldehyde [%]
    "RPM": (100, 300)     # Agitation speed
}

# 2) How many factors?
k = len(factors)

# 3) Choose CCD type
#    Face-centered (CCF): star points on cube faces, practical with hard bounds
#    Circumscribed (CCC): star points outside cube, rotatable option available
face_type = 'ccf'  # face-centered default
alpha_type = 'o'   # orthogonal (alpha parameter ignored for ccf)

# For rotatable circumscribed design, uncomment:
# face_type = 'ccc'
# alpha_type = 'r'

# 4) Center point replicates (factorial block, star block)
#    Total center points = factorial + star
center = (2, 3)  # 2 in factorial + 3 in star = 5 total center points

# 5) Build the CCD in coded units
#    Columns will match the insertion order of the keys in `factors`.
ccd_coded = ccdesign(k, center=center, face=face_type, alpha=alpha_type)
ccd_coded_df = pd.DataFrame(ccd_coded, columns=list(factors.keys()))

Note: In pyDOE3, the center=(center1, center2) parameter assigns center points to two separate blocks: center1 goes into the factorial block and center2 into the axial (star) block. This design choice reflects real-world experimental workflows. Typically, you start with a factorial design and include one or two center points to check for curvature. If curvature is detected, you extend the design to a full CCD by adding axial runs and additional center points (often 3–4 more) to improve the estimate of pure error and stabilize the response surface model. By separating the center point counts, pyDOE3 lets you mimic this sequential experimentation strategy.

What you get:

• Factorial points at ±1
• Center points at 0
• Star points at ±alpha on one factor with all others at 0

Run	T_C	P_bar	CoF_pct	RPM
0	-1.0	-1.0	-1.0	-1.0
1	-1.0	-1.0	-1.0	1.0
2	-1.0	-1.0	1.0	-1.0
3	-1.0	-1.0	1.0	1.0
4	-1.0	1.0	-1.0	-1.0
5	-1.0	1.0	-1.0	1.0
6	-1.0	1.0	1.0	-1.0
7	-1.0	1.0	1.0	1.0
8	1.0	-1.0	-1.0	-1.0
9	1.0	-1.0	-1.0	1.0
10	1.0	-1.0	1.0	-1.0
11	1.0	-1.0	1.0	1.0
12	1.0	1.0	-1.0	-1.0
13	1.0	1.0	-1.0	1.0
14	1.0	1.0	1.0	-1.0
15	1.0	1.0	1.0	1.0
16	0.0	0.0	0.0	0.0
17	0.0	0.0	0.0	0.0
18	-1.0	0.0	0.0	0.0
19	1.0	0.0	0.0	0.0
20	0.0	-1.0	0.0	0.0
21	0.0	1.0	0.0	0.0
22	0.0	0.0	-1.0	0.0
23	0.0	0.0	1.0	0.0
24	0.0	0.0	0.0	-1.0
25	0.0	0.0	0.0	1.0
26	0.0	0.0	0.0	0.0
27	0.0	0.0	0.0	0.0
28	0.0	0.0	0.0	0.0

Step 3: Map coded values to real-world settings

Mapping coded units to real units keeps the plan practical and readable. The linear mapping below works for both face-centered and rotatable CCDs.

Mathematically, for each factor with lower bound L and upper bound H, coded values x_coded in [−1, 1] map to real units via a simple linear transformation:

x_real = (L + H)/2 + x_coded × (H − L)/2

Now apply that mapping in Python and add helpful labels to each point:

def coded_to_real(df_coded: pd.DataFrame, spec: dict) -> pd.DataFrame:
    df = df_coded.copy()
    for name, (low, high) in spec.items():
        mid = (low + high) / 2.0
        half_range = (high - low) / 2.0
        df[name] = mid + df[name] * half_range
    return df

ccd_real = coded_to_real(ccd_coded_df, factors)

# Optional: classify point type for clarity
def classify_point(row) -> str:
    vals = np.abs(row.values)
    if np.allclose(vals, 0):
        return "Center"
    # For face-centered, axial points are at ±1 on one factor, 0 on others
    # For circumscribed, they're at ±α (where α > 1 typically)
    is_axial = (vals.max() > 1.0) and (np.sum(vals > 0.1) == 1)  # One non-zero factor
    if is_axial:
        return "Axial"
    is_factorial = np.all(np.isin(np.round(vals, 6), [1.0]))
    if is_factorial:
        return "Factorial"
    return "Other"

point_types = ccd_coded_df.apply(classify_point, axis=1)
ccd_real.insert(0, "PointType", point_types)

print(ccd_real.head())

Example output (first rows) will include factorial, axial, and center labeled points with real-world values:

Run	PointType	T_C	P_bar	CoF_pct	RPM
0	Factorial	20.0	1.0	2.0	100.0
1	Factorial	20.0	1.0	2.0	300.0
2	Factorial	20.0	1.0	6.0	100.0
3	Factorial	20.0	1.0	6.0	300.0
4	Factorial	20.0	3.0	2.0	100.0

Step 4: Randomize and finalize the plan

Randomization protects against time trends and other lurking variables. Center point replicates are already built into the design; you can also add full design replicates if measurement error is high.

# Randomize run order for execution
randomized = ccd_real.sample(frac=1, random_state=42).reset_index(drop=True)

# Add a readable Run index
randomized.index += 1
randomized.index.name = "Run"

print(randomized.head())

Output:

Run	PointType	T_C	P_bar	CoF_pct	RPM
1	Center	30.0	2.0	4.0	200.0
2	Center	30.0	2.0	4.0	200.0
3	Factorial	40.0	3.0	2.0	100.0
4	Other	30.0	2.0	2.0	200.0
5	Factorial	40.0	1.0	2.0	100.0

Step 5: Export to Excel

Export a clean, randomized plan for execution or review.

randomized.to_excel("ccd_plan.xlsx")  # index label "Run" is included

That’s it! Have fun ;D

Up next:

<< Create a Fractional Factorial Design in Python >>

<< Understanding Bayesian Optimization >>