本文是《Quality by Design for ANDAs: An Example for Immediate-Release Dosage Forms》第一個處方的R語言解決方案。

第一個處方研究評估原料藥粒徑分布、MCC/Lactose比例、崩解劑用量對制劑CQAs的影響。

第二處方研究用于理解顆粒外加硬脂酸鎂和滑石粉對片劑質量和可生產性的影響。處方研究在實驗規模進行（1.0kg，5000片）。

處方研究1的目的是選擇MCC/Lactose的比例以及崩解劑的用量，并理解這些變量是否與原粒藥粒徑分布有交互作用。這個研究也用來確立處方穩健性。這些處方因素對響應變量的影響的帶3個中心點的23全因子DOE列于表21。

崩解劑加到顆粒內，用量范圍為1-5%。MCC/Lactose的比例為1：2，1：1，2：1。原料藥的重量占10%。滑石粉的用量占2.5%。硬脂酸鎂的用量占1%。片重為200.0mg。

用三種滾壓力滾壓物料。調整壓片力（允許的范圍是11-13.0kP）以研究片劑硬度為12.0kP時的溶出。選擇12.0kP的片劑硬度以研究處方變量對溶出的影響因為更高的硬度是溶出的最差條件。如果固定壓片力則可以混雜片劑硬度的影響。

表21列出了考察的因子和響應。

表22列出了溶出，含量均一性，粉末流動系數和10KN壓片力的片劑硬度。

library(FrF2)

study1<-FrF2(nruns=8,nfactors=3, ncenter=3,factor.names=list(A=c("-1","+1"),B=c("-1","+1"),C=c("-1","+1")), replications=1,randomize=FALSE)

> study1

y1<-c(99.5,77.0,98.7,86.0,99.0,76.0,99.0,84.0,91.0,89.4,92.0)

y3<-c(4.1, 2.9, 4.0, 3.2, 5.0, 3.8, 5.1, 4.0, 4.0, 3.9, 4.1)

y5<-c(6.16,8.46,6.09,8.46,4.97,7.56,4.77,7.25,6.62,6.66,6.46)

y7<-c(9.1,8.3,? 9.1, 8.6, 13.5,12.5,12.9,13.2,10.6,10.9,11.3)

study1 <- add.response( study1, y1)

study1 <- add.response( study1, y3)

study1 <- add.response( study1, y5)

study1 <- add.response( study1, y7)

mod1 <- lm( y1 ~A*B, data = study1)

anova(mod1)

將Residuals分解為彎曲性、純誤差、Lack of Fit:

SS Pure Error=((91.0)^2+(92.0)^2+(89.4)^2-(91.0+92.0+89.4)^2/3)=3.44

純誤差的自由度為3-1=2

MS Pure Error =3.44/2=1.72

SS Curvature=8*3*((84.0+98.7+77.0+99.0+99.0+99.5+86.0+76.0)/8-(91.0+92.0+89.4)/3)^2/(8+3)= 1.77

彎曲性的自由度為1

MS Curvature=1.77/1=1.77

F Curvature=MS Curvature/MS Residuals = 1.74

SS Lack of fit=7.88-3.44-1.77=2.67

Lack of fit的自由度為4

MS lack of fit= SS Lack of fit/4=0.67

F lack of fit=MS lack of fit/MS Pure Error=0.39

Pr(>F) Curvature=1-pf(1.74,1,6)=0.2352431

Pr(>F) lack of fit=1-pf(0.39,4,2)=0.8079788

最終的方差分析表如下：

Response: y1

???????? ??????Df Sum Sq Mean Sq F value??? Pr(>F)???

A???????? ??????1 669.78? 669.78 ?657.72 ?2.316471e-07 ***

B???????? ??????1? 32.81?? 32.81 ??32.22 ?0.001287723 **

A:B????? ???????1? 39.60?? 39.60 ??38.89 ?0.00078691 ?***

Curvature????? ??1?? 1.77??? 1.77??? 1.74? 0.2352431

Residuals????? ??6?? 6.11??? 1.02???? -????? -

lack of fit? ??????4?? 2.67??? 0.67??? 0.39? 0.8079788?

Pure Error ???????2?? 3.44??? 1.72???

表23是溶出度的ANOVA分析結果。

如表23所示，溶出度的彎曲效應并不顯著。用所有的數據（包括中心點）擬合模型。

library(FrF2)

study1<-FrF2(nruns=8,nfactors=3, ncenter=3,factor.names=list(A=c("-1","+1"),B=c("-1","+1"),C=c("-1","+1")), replications=1,randomize=FALSE)

study1

y1<-c(99.5,77.0,98.7,86.0,99.0,76.0,99.0,84.0,91.0,89.4,92.0)

y3<-c(4.1, 2.9, 4.0, 3.2, 5.0, 3.8, 5.1, 4.0, 4.0, 3.9, 4.1)

y5<-c(6.16,8.46,6.09,8.46,4.97,7.56,4.77,7.25,6.62,6.66,6.46)

y7<-c(9.1,8.3,? 9.1, 8.6, 13.5,12.5,12.9,13.2,10.6,10.9,11.3)

study1 <- add.response( study1, y1)

study1 <- add.response( study1, y3)

study1 <- add.response( study1, y5)

study1 <- add.response( study1, y7)

mod1 <- lm( y1 ~A*B, data = study1)

anova(mod1)

> anova(mod1)

Analysis of Variance Table

Response: y1

????????? Df Sum Sq Mean Sq F value??? Pr(>F)???

A????????? 1 669.78? 669.78 595.188 4.951e-08 ***

B????????? 1? 32.81?? 32.81? 29.152 0.0010093 **

A:B??????? 1? 39.60?? 39.60? 35.194 0.0005802 ***

Residuals? 7?? 7.88??? 1.13?????????????????????

---

Signif. codes:? 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

將Residuals分解為純誤差、Lack of Fit:

SS Pure Error=((91.0)^2+(92.0)^2+(89.4)^2-(91.0+92.0+89.4)^2/3)=3.44

純誤差的自由度為3-1=2

MS Pure Error =3.44/2=1.72

SS Lack of fit=7.88-3.44=4.44

Lack of fit的自由度為7-2=5

MS lack of fit= SS Lack of fit/5=0.888

F lack of fit=MS lack of fit/MS Pure Error=0.516

Pr(>F) lack of fit=1-pf(0.516,5,2)=0.7618

library(daewr)

fullnormal(coef(mod1)[-1], alpha=.025)

effects <-coef(mod1)

effects <-effects[2:4]

effects <-effects[ !is.na(effects) ]

library(daewr)

halfnorm(effects, names(effects), alpha=.25)

library(BsMD)

LenthPlot(mod1, main = "Lenth Plot of Effects")

A.num <-study1$A

levels(A.num) <- c(10,30)

B.num <- study1$B

levels(B.num) <- c(1,5)

A.num <- as.numeric(as.character(A.num))

B.num <- as.numeric(as.character(B.num))

mod1 <- lm( y1 ~A.num*B.num, data = study1)

library(rsm)

contour(mod1, ~ A.num+B.num)

persp(mod1, ~ A.num+B.num, zlab="y1", contours=list(z="bottom"))

如圖10的半正態圖和表24的未調整模型的ANOVA結果所示，影響溶出度的顯著因子為A（原料藥粒徑）,B（崩解劑用量），AB（原料粒徑與崩解劑用量交互作用）。

圖11展示原粒藥粒徑與崩解劑用量對溶出度的影響。

A <- factor(c(1,-1, 1,-1,-1,-1,1, 1))

B <- factor(c(1, 1,-1,-1, 1,-1,1,-1))

C <- factor(c(1,-1,-1, 1, 1,-1,1, 1))

y1<-c(84.0, 98.7, 77.0, 99.0, 99.0, 99.5, 86.0,76.0)

y3<-c(4.0,? 4.0,? 2.9,? 5.0,? 5.1,? 4.1,? 3.2, 3.8 )

y5<-c(7.25, 6.09, 8.46, 4.97, 4.77, 6.16, 8.46,7.56)

y7<-c(13.2, 9.1,? 8.3,? 13.5, 8.3,? 9.1,? 8.6, 12.5)

study1<- data.frame (A=A,B=B,C=C,y1=y1, y3=y3, y5=y5, y7=y7)

write.csv( plan, file = "E://data//study1.csv", row.names = FALSE )

study1 <- read.csv("E://data//study1.csv")

mod1 <- lm( y1 ~A*B, data = study1)

anova(mod1)

> anova(mod1)

Analysis of Variance Table

Response: y1

????????? Df Sum Sq Mean Sq? F value?? Pr(>F)???

A????????? 1 669.78? 669.78 1003.416 5.92e-06 ***

B????????? 1? 32.80?? 32.80?? 49.146 0.002180 **

A:B??????? 1? 39.60?? 39.60?? 59.333 0.001529 **

Residuals? 4?? 2.67??? 0.67?????????????????????

---

Signif. codes:? 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

with(study1, (interaction.plot(A,B,y1, type = "b", pch = c(18,24,22), leg.bty = "o", main = "Interaction Plot of A and B", xlab = "A",ylab = "y1")))

#以下將因子變量轉為自然變量

A.num <-study1$A

levels(A.num) <- c(10,30)

B.num <- study1$B

levels(B.num) <- c(1,5)

A.num <- as.numeric(as.character(A.num))

B.num <- as.numeric(as.character(B.num))

mod1 <- lm( y1 ~A.num*B.num, data = study1)

contour(mod1, ~ A.num + B.num)

persp(mod1, ~ A.num + B.num, zlab=" y1", contours=list(z="bottom"))

影響崩解時間的顯著因子

崩解劑是唯一影響片劑崩解時間的顯著因子。但是所有的批次都在小于4分鐘內快速的崩解。

影響片劑含量的顯著因子

所有的批次含量均在標準范圍內(95.0-105.0% w/w) ，沒有因子對含量有顯著影響。