クリアランスとメカニズム

ファーマコキネティクス―演習による理解
杉山 雄一,山下 伸二,加藤 基浩
南山堂 / ￥ 6,300 ()
在庫あり。

設問1

> data <- read.csv("/Users/kzfm/PK/pk33.csv")
> data
  drag  dose      AUC   fe
a    A 1e+00    169.0 0.13
b    A 1e+01   1720.0 0.14
c    A 1e+02  22400.0 0.18
d    A 1e+03 530000.0 0.42
e    B 1e-02      2.4 0.17
f    B 1e-01     23.4 0.17
g    B 1e+00    167.0 0.19
h    B 1e+01    658.0 0.25

> attach(data)
> CLtot <- dose/AUC * 1000
> CLtot
[1]  5.917160  5.813953  4.464286  1.886792  4.166667  4.273504  5.988024 15.197568

> CLr <- CLtot * fe
> CLr
[1] 0.7692308 0.8139535 0.8035714 0.7924528 0.7083333 0.7264957 1.1377246 3.7993921

> CLh <- CLtot - CLr
> CLh
[1]  5.147929  5.000000  3.660714  1.094340  3.458333  3.547009  4.850299 11.398176

設問2

CLr = fb * GFR つまりfb = VLr/GFR

> GFR <- 8
> QH <- 60

> fb_a <- CLr[1]/8
> fb_b <- CLr[5]/8
> fb_a
[1] 0.09615385
> fb_b
[1] 0.08854167

CLh = (QH*fb*CLh_uint)/(QH+fb*CLh_uint)つまりCLh_uint=CLh*QH/(fb*(QH-CLh))

> CLh_uint_a <- CLh[1]*QH/(fb_a*(QH-CLh[1]))
> CLh_uint_a
[1] 58.56311
> CLh_uint_b <- CLh[5]*QH/(fb_b*(QH-CLh[5]))
> CLh_uint_b
[1] 41.44783

設問3

DrugA

CLrの投与量変化が見られないことからfbの変化はない。つまり投与量増加に伴うCLhの減少は CLh_uintの低下、つまり代謝過程の飽和

DrugB

投与量増加に伴い、CLh,CLrが増加。CLrの増加はfbの増加によるもので、これはつまり血中蛋白結合の 飽和をあらわしている。 肝固有クリアランスとの関係はなんとなく。

繰り返し経口投与の続き

ファーマコキネティクス―演習による理解
杉山 雄一,山下 伸二,加藤 基浩
南山堂 / ￥ 6,300 ()
在庫あり。

設問 4

肝クリアランスと肝血流の大きさから律速を判断。

v = Vmax*[S] / (Km + [S])

肝クリアランスに飽和が見られるから

CLh_int = Vmax/Km * 1/Cu

ここでCu = Cpss * fp

vh = Cb*CLh vr = Cpss*Rb*CLr

設問 5

v=D/t*(F-fe),C=Cpss*fp でMicaehlis-Menten式の形でフィッティング

> data <- read.csv("/Users/kzfm/PK/pk32b.csv")
> data
  Dose      C     v
a   10  0.212 0.308
b   30  0.675 0.850
c  100  2.400 2.420
d  300  7.850 5.750
e  500 13.800 7.920
> attach(data)
> e <- nls(v ~ Vmax*C/(Km+C),start=list(Vmax=10,Km=10))
> e
Nonlinear regression model
  model:  v ~ Vmax * C/(Km + C) 
   data:  parent.frame() 
 Vmax    Km 
15.09 12.59 
 residual sum-of-squares: 0.01301

繰り返し経口投与

ファーマコキネティクス―演習による理解
杉山 雄一,山下 伸二,加藤 基浩
南山堂 / ￥ 6,300 ()
在庫あり。

設問1

省略

設問 2

CL = (F * (Dose/t)) / Cp,ss 10mgのデータからFを決定する

CLb=CLp/RbからRbをつかって全血基準のクリアランスを求める

t1/2をもとめる

一般に、薬物を持続投与する場合、消失半減期の5倍の時間が経っている時には定常状態に達しているとみてよいので これをもって定常状態の判断に利用してよい

設問 3

> data <- read.csv("/Users/kzfm/PK/pk3-2.csv")
> data
  dose     cpss   fe
A   10 0.000424 0.63
B   30 0.001350 0.66
C  100 0.004800 0.71
D  300 0.015700 0.77
E  500 0.027500 0.81

> attach(data)
> fp <- 0.5
> Rb <- 0.7
> Vd <- 30
> ka <- 1.5

> CLb <- (dose/12)/(cpss*1000*Rb)
> CLb
[1] 2.807727 2.645503 2.480159 2.274795 2.164502

> CLr <- CLb*fe
> CLr
[1] 1.768868 1.746032 1.760913 1.751592 1.753247

> CLh <- CLb - CLr
> CLh
[1] 1.0388589 0.8994709 0.7192460 0.5232029 0.4112554

設問1はクリアランスのパラメータを求める

設問2,3はi.v., p.o.それぞれの状況でパラメータを変えた場合どう変化するかを推定する

ファーマコキネティクス―演習による理解
杉山 雄一,山下 伸二,加藤 基浩
南山堂 / ￥ 6,300 ()
在庫あり。

設問1

> data <- read.csv("/Users/kzfm/PK/pk31.csv")
> data
  CLtot   Vd   fe  fb  Ka
A   0.2  0.1 1.00 0.1 1.5
B   0.2 10.0 1.00 0.1 1.5
C   1.3  0.1 0.02 0.1 1.5
D   1.3 10.0 0.02 0.1 1.5

> Qh <- 1.4
> Qr <- 1.1
> GFR <- 0.17
> Vb <- 0.077

> CLh <- (1-fe)*CLtot
> CLh
[1] 0.000 0.000 1.274 1.274

> CLr <- fe*CLtot
> CLr
[1] 0.200 0.200 0.026 0.026

> CLh_uint <- Qh*CLh/(fb*(Qh-CLh))
> CLh_uint
[1]   0.0000   0.0000 141.5556 141.5556

> CLr_uint <- Qr*(CLr -fb*GFR)/((Qr-CLr+fb*GFR)*fb)
> CLr_uint
[1] 2.19520174 2.19520174 0.09074244 0.09074244

設問2

AUCiv = D/CLtot

1. fbの変動によって影響を受けるパラメータはVd,CLtot
2. Qhの変動により影響を受けるパラメータはCLh
   • CLh = Qh*fb*CLh_uint/(Qh + fb*CLh_uint)なので分母の項の大小で律速が決まる
3. CLr_uintの変動により影響をうけるのはCLr
4. CLh_uintの変動により影響をうけるのはCLh

設問3

flip-flopを考える

Tweeting from R

> library("RCurl")
> opts <- curlOptions(header = FALSE, 
+   userpwd = "[user]:[pass]", netrc = FALSE)
> tweet <- function(status){
+   method <- "http://twitter.com/statuses/update.xml?status="
+   encoded_status <- URLencode(status)
+   request <- paste(method,encoded_status,sep = "")
+   postForm(request,.opts = opts)
+ }
> xx<-rnorm(2)
> tweet(paste(xx,collapse="|"))

つぶやいた

SVMかなんかを利用して分類結果をつぶやくボットをつくっても面白いかも。

モーメント解析。AUCとAUMCをもとめるのにPKというライブラリを利用した。

ファーマコキネティクス―演習による理解
杉山 雄一,山下 伸二,加藤 基浩
南山堂 / ￥ 6,300 ()
通常2～4週間以内に発送

設問1

library(PK)
time <- c(0.25, 0.5, 1, 2, 4, 6, 8, 10, 12, 16)
l_po_sol <- c(25.92, 33.75, 42.08, 28.81, 11.26, 4.24, 1.76, 0.59, 0.22)
l_po_tab <- c(1.3, 4.2, 11.33, 20.35, 20.15, 13.17, 6.66, 4.57, 1.81, 0.67)
l_iv     <- c(91.21, 76.53, 61.46, 32.36, 11.9, 4.67, 2.15, 0.67, 0.25)
u_po_sol <- c(25.95, 33.75, 42.08, 28.81, 11.26)
u_po_tab <- c(11.33, 20.35, 20.15, 13.17)
u_iv     <- c(91.21, 76.53, 61.46, 32.36, 11.9)

auc(time[1:9],l_po_sol,n.tail=3)
Estimation for a serial sampling design

                  AUC     AUMC
observed     113.8863 1209.474
interpolated       NA       NA
infinity     140.1590 3005.290

auc(time,l_po_tab,n.tail=3)
Estimation for a serial sampling design

                  AUC      AUMC
observed     106.8875  953.0328
interpolated       NA        NA
infinity     125.3667 1499.8264

auc(time[1:9],l_iv,n.tail=3)
Estimation for a serial sampling design

                  AUC     AUMC
observed     193.5675 4534.682
interpolated       NA       NA
infinity     198.9328 5139.201

クリアランスは投与量をAUCでわればよい。

#CLtot = DOSE / AUC

0.1 / 198.9328 * 1000000 # = 502.6823

#Fsol = AUCpo(sol) / AUCiv
140.1590 / 198.9328 # 0.7045545

#Ftab = AUCpo(tab) / AUCiv
125.3667 / 198.9328 # 0.6301962

MRTはAUMCをAUCでわればよいはずなのだが、何かおかしい。演習の答えとあわん。

#MRT = AUMC / AUC
#Vss = CL * MRT

設問2

デコンボリューションはまたあとで。

09.06.02

超ぼけてた

昨晩AUMCがおかしかったのは、timeとconcの順番が逆だったという凡ミス。auc.completeでやれば明示的に指定できるのでこっちで書いてみた。

library(PK)
time <- c(0.25, 0.5, 1, 2, 4, 6, 8, 10, 12, 16)
l_po_sol <- c(25.92, 33.75, 42.08, 28.81, 11.26, 4.24, 1.76, 0.59, 0.22)
L_po_tab <- c(1.3, 4.2, 11.33, 20.35, 20.15, 13.17, 6.66, 4.57, 1.81, 0.67)
l_iv     <- c(91.21, 76.53, 61.46, 32.36, 11.9, 4.67, 2.15, 0.67, 0.25)

l_po_sol_parm <- auc.complete(conc=c(0,l_po_sol), time=c(0,time[1:9]), n.tail=3)
l_po_tab_parm <- auc.complete(conc=c(0,l_po_tab), time=c(0,time), n.tail=3)
l_iv_parm <- auc.complete(conc=c(l_iv), time=c(time[1:9]), n.tail=3)

MRTをもとめる(AUMC / AUC)

> l_iv_parm$est[3,2]/l_iv_parm$est[3,1]
[1] 2.117064
> l_po_sol_parm$est[3,2]/l_po_sol_parm$est[3,1]
[1] 2.421341
> l_po_tab_parm$est[3,2]/l_po_tab_parm$est[3,1]
[1] 5.151322

これでうまくいきそう。

ファーマコキネティクス―演習による理解
杉山 雄一,山下 伸二,加藤 基浩
南山堂 / ￥ 6,300 ()
通常2～4週間以内に発送

設問1

フィッティングのための初期値を求める。

t <- c(2.1, 4.1, 6.1, 8.1, 10.1, 14.1, 18.1, 22.1, 30.1, 45.2, 60.2, 
120.3, 180.5, 241.8, 470, 948, 1423, 2187, 2831, 4305)
conc <- c(20.5, 17.5, 14.5, 12.5, 13.0, 12.0, 11.0, 9.1, 9.6, 5.6, 4.9,
3.2, 2.0, 1.8, 0.90, 0.85, 0.70, 0.45, 0.43, 0.39)

logc <- log(conc)

lm(logc[15:20] ~ t[15:20])

Call:
lm(formula = logc[15:20] ~ t[15:20])

Coefficients:
(Intercept)     t[15:20]  
 -0.0399743   -0.0002441  

exp(-0.0399)
[1] 0.9608855

beta <- function(t){exp(-0.0399743-t)}
ac <- conc[1:14] - beta(conc[1:14])

lm(ac ~ t[1:14])

Call:
lm(formula = ac ~ t[1:14])

Coefficients:
(Intercept)      t[1:14]  
   13.19420     -0.06273  

設問2

2,3exp式でフィッティングし、あてはめの良さをくらべる。

e2 <- nls(conc ~ a*exp(-alpha*t)+b*exp(-beta*t),
start=list(a=16,alpha=0.028,b=2,beta=0.000637))

e3 <- nls(conc ~ a*exp(-alpha*t)+b*exp(-beta*t)+c*exp(-gamma*t),
start=list(a=16,alpha=0.13,b=8.4,beta=0.0114,c=0.96,gamma=0.000239))

あてはめのよさは赤池の情報量基準で。Rの場合はAIC関数でよい。

AIC(e2) # [1] 56.79388
AIC(e3) # [1] 44.93942

というわけで、exp3式のほうがあてはまりがよい。

設問3

skip

ファーマコキネティクス―演習による理解
杉山 雄一,山下 伸二,加藤 基浩
南山堂 / ￥ 6,300 ()
在庫あり。

微分方程式を解くためにodesolveを利用する。初期値の与えかたとか、ぐるぐる回していくイメージをつかむのに微妙に時間がかかったが、覚えてしまうと簡単な方程式を解かせたりシミュレーションするのには楽になる。

設問1

吸収過程を考えないコンパートメントモデル

library(odesolve)
params <- c(kel = 0.037)
times <- c(0,(1:40))
dydt <- function(t,y,p){
  kel <- p['kel']
  y <- -kel*y
  list(c(y))
}
y <- lsoda(c(y = 8.6),times,dydt,params)
plot(y,ylab="Conc",xlab="time (hour)")

設問2&3

今度は吸収過程を考える。設問3は吸収速度定数を1/10にした時にどういう挙動を示すか。

params <- c(kel = 0.037,ka = 0.1, D = 5000, V = 580)
times <- c(0,(1:120))
dydt <- function(t,y,p){
  kel <- p['kel']
  ka <- p['ka']
  D <- p['D']
  V <- p['V']
  y <- ka * (D/V) * exp(-ka*t) -kel*y 
  list(c(y))
}
y <- lsoda(c(y = 0),times,dydt,params)

params2 <- c(kel = 0.037,ka = 0.01, D = 5000, V = 580)
y2 <- lsoda(c(y = 0),times,dydt,params2)

plot(y,ylim=c(0,6),ylab="Conc",xlab="time(hour)")
par(new=T)
plot(y2,ylim=c(0,6),axes=F,xlab="",ylab="")

今回バイオアバイラビリティ(F)を1とおいているが、これを変化させた時の挙動も気になるのでやってみた。設問の式にFを組み込んでF=0.1とした時の結果

dydt2 <- function(t,y,p){
  kel <- p['kel']
  ka <- p['ka']
  D <- p['D']
  V <- p['V']
  F <- p['F']
  y <- ka * (D/V) * F * exp(-ka*t) -kel*y 
  list(c(y))
}

params3 <- c(kel = 0.037,ka = 0.1, D = 5000, V = 580, F = 0.1)
y3 <- lsoda(c(y = 0),times,dydt2,params3)

plot(y2,ylim=c(0,3),axes=F,xlab="time(hour)",ylab="Conc")
par(new=T)
plot(y3,ylim=c(0,3),axes=F,xlab="",ylab="")

あーなるほど、そうだよなぁと一人で妙に納得したのであった。

ファーマコキネティクス―演習による理解
杉山 雄一,山下 伸二,加藤 基浩
南山堂 / ￥ 6,300 ()
在庫あり。

odesolveを覚えるためにロトカ＝ヴォルテラの方程式サンプルを。

lvmodel <- function(t, x, parms) {
s <- x[1] # substrate
p <- x[2] # producer
k <- x[3] # consumer
with(as.list(parms),{
import <- approx(signal$times, signal$import, t)$y
ds <- import - b*s*p + g*k
dp <- c*s*p - d*k*p
dk <- e*p*k - f*k
res<-c(ds, dp, dk)
list(res)
})
}
## vector of timesteps
times <- seq(0, 100, length=101)
## external signal with rectangle impulse
signal <- as.data.frame(list(times = times,
import = rep(0,length(times))))
signal$import[signal$times >= 10 & signal$times <=11] <- 0.2
## Parameters for steady state conditions
parms <- c(a=0.0, b=0.0, c=0.1, d=0.1, e=0.1, f=0.1, g=0.0)
## Start values for steady state
y<-xstart <- c(s=1, p=1, k=1)
## Classical RK4 with fixed time step
out1 <- as.data.frame(rk4(xstart, times, lvmodel, parms))
## LSODA (default step size)
out2 <- as.data.frame(lsoda(xstart, times, lvmodel, parms))
## LSODA: with fixed maximum time step
out3 <- as.data.frame(lsoda(xstart, times, lvmodel, parms, hmax=1))
par(mfrow=c(2,2))
plot (out1$time, out1$s, type="l", ylim=c(0,3))
lines(out2$time, out2$s, col="red", lty="dotted")
lines(out3$time, out3$s, col="green", lty="dotted")
plot (out1$time, out1$p, type="l", ylim=c(0,3))
lines(out2$time, out2$p, col="red", lty="dotted")
lines(out3$time, out3$p, col="green", lty="dotted")
plot (out1$time, out1$k, type="l", ylim=c(0,3))
lines(out2$time, out2$k, col="red", lty="dotted")
lines(out3$time, out3$k, col="green", lty="dotted")
plot (out1$p, out1$k, type="l")
lines(out2$p, out2$k, col="red", lty="dotted")
lines(out3$p, out3$k, col="green", lty="dotted")

参考

• R で微分方程式：odesolve
• Rで常微分方程式を解いてみた．Michaelis-Menten kineticsを例として
• Rの練習のために力学系を描いてみた