_ | 覦覈襦 | 豕蠏手 | 殊螳 | 譯殊碁
FrontPage › 觜蠏覿

Contents

1
2 n谿 curve fitting
3 self-starting function
4 inverse.predict
5 bi-exponential
6 exp * 3
7 log-normal function
8 pareto
9 覃燕 3螳
10 譯手鍵 願鍵
11 smooth
12 觜 蠏 襤郁規螳
13 beta
14 谿瑚襭



library(nlstools)
plotfit(fit, smooth=TRUE)
overview(fit)

fit.boot <- nlsBoot(fit, niter = 200)
summary(fit.boot)

1 #

れ 襭襯 覲企, .
install.packages("car")
library("car")
plot(population ~ year, data=USPop, main="(a)")
abline(lm(population ~ year, data=USPop))
nls01.png

2 n谿 curve fitting #

curve fitting 企慨.
xx <- USPop$year
fit <- lm(population ~ poly(year, 2), data=USPop)
lines(xx, predict(fit, data.frame(year=xx)), col="purple")
nls02.png

蠏碁磯 USPop 一危磯 螳 logistic growth model 企手 .
curve(1/(1+exp(-(-1 + 1*x))), from=-5, to=5, main="(b)")
abline(h=1/2, lty=2)
abline(v=1, lty=2)
nls03.png

3 self-starting function #

n谿襦 fitting 蟆 螳 覃 螻, 襷 覿覿 れ 螻旧れ 蟆企. self-starting functionれ企.

nls04.png

n谿 curve fitting企 螳 self-starting function 伎覃 蟆企.
fit <- nls(population ~ SSlogis(year, Asym, xmid, scal), data=USPop)
summary(fit)
xx <- USPop$year
plot(population ~ year, data=USPop, main="(a)")
lines(xx, predict(fit, data.frame(year=xx)), col="purple")
nls05.png

SSasymp(input, Asym, R0, lrc)
SSasympOff(input, Asym, lrc, c0)
SSasympOrig(input, Asym, lrc)
SSbiexp(input, A1, lrc1, A2, lrc2)
SSfol(Dose, input, lKe, lKa, lCl)
SSfpl(input, A, B, xmid, scal)
SSgompertz(input, Asym, b2, b3)
SSlogis(input, Asym, xmid, scal)
SSmicmen(input, Vm, K)
SSweibull(input, Asym, Drop, lrc, pwr)

self-starting function 覈 螳襯 蟆 谿碁 觜蟲企蓋 蟆企. self-starting function 朱誤磯 麹(襯 れ, help 螳) 覲襯 j, 螳 覃 蟇郁, 蟇郁.
xx <- USPop$year
plot(population ~ year, data=USPop, main="(a)")

fit1 <- nls(population ~ SSlogis(year, Asym, xmid, scal), data=USPop)
lines(xx, predict(fit1, data.frame(year=xx)), col="purple")

fit2 <- nls(population ~ SSweibull(year, Asym, Drop, lrc, pwr), data=USPop)
lines(xx, predict(fit2, data.frame(year=xx)), col="red")

fit3 <- nls(population ~ SSfpl(year, A, B, xmid, scal), data=USPop)
lines(xx, predict(fit3, data.frame(year=xx)), col="green")
nls06.png

4 inverse.predict #

install.packages("chemCal")
library("chemCal")

fit1 <- rlm(dau ~ x, data=tmp)
inverse.predict(fit1, 500000) #dau螳 500000 る x 覈伎伎 螳?

5 bi-exponential #

> head(df)
  seq    u
1   1 4311
2   2 2381
3   3 1981
4   4 1529
5   5 1386
6   6 1340

x1 <- df[1:7,]
x2 <- df[8:nrow(df),]

cc1 <- coef(lm(log(u) ~ seq, data=x1))
cc2 <- coef(lm(log(u) ~ seq, data=x2))

fit <- nlsLM(u ~ a*exp(b * seq) + c * exp(d * seq),
             start = list(a=cc1[1], b=cc1[2], c=cc2[1], d=cc2[2]), 
             data=df)

6 exp * 3 #

#exp * 3
r <- sort(unique(df.ko$ranking))
xlim <- c(1,max(r))
ylim <- c(0, max(df.ko$amt) + max(df.ko$amt) * 0.3) 
plot(df.ko$ranking, df.ko$amt, cex=0.1, pch=19, xlim=xlim, ylim=ylim, xlab="谿", ylab="u", main="bi-exp")

x1 <- df.ko[df.ko$ranking <= 7,]
x2 <- df.ko[df.ko$ranking >  7 & df.ko$ranking <= 20 ,]
x3 <- df.ko[df.ko$ranking >  21,]

cc1 <- coef(lm(log(amt) ~ ranking, data=x1))
cc2 <- coef(lm(log(amt) ~ ranking, data=x2))
cc3 <- coef(lm(log(amt) ~ ranking, data=x3))

fit <- nlsLM(u ~ a*exp(b * seq) + c*exp(d * seq) + e*exp(f * seq),
             start = list(a=cc1[1], b=cc1[2], c=cc2[1], d=cc2[2], e=cc3[1], f=cc3[2]), 
             data=df)
summary(fit)
lines(r, predict(fit, data.frame(ranking=r)), col="red")
#text(df$seq, df$u, df$w, cex=0.5, pos=4)

amt1 <- predict(fit, data.frame(ranking=r))
amt1 / sum(amt1)

r <- sort(unique(df.ko$ranking))
xlim <- c(1,max(r))
ylim <- c(0, max(df.ko$amt) + max(df.ko$amt) * 0.3) 
plot(df.ko$ranking, df.ko$amt, cex=0.1, pch=19, xlim=xlim, ylim=ylim, xlab="", ylab="襷れ", main="蟲")
lines(r, predict(fit, data.frame(ranking=r)), col="red")

lines(r, coef(fit)[1]*exp(coef(fit)[2] * df$seq), col="blue")
lines(r, coef(fit)[3]*exp(coef(fit)[4] * df$seq), col="green")
lines(r, coef(fit)[5]*exp(coef(fit)[6] * df$seq), col="black")



7 log-normal function #

fun <- function(x, mu, sd){
    a <- 1/(x * sqrt(2*pi*sd))
    b <- (log(x)-mu)^2 * -1
    c <- 2 * sd^2
    return (a * exp(b/c))
}


#data
val <- c(462604.2114,724545.9975,677961.1773,556621.1925,572068.4823,550408.0245,442748.9577,463050.285,485241.4509,406399.9335,394768.1661,474606.3792,408121.4988,254290.8273,234755.1933,201690.9834,
175243.2,171106.0665,228024.2613,199481.5251,141505.8969,148199.988,138756.7692,140078.0631,146512.2765,183688.7274,169084.7955,103705.1421,107350.3998,106554.8355,99301.161,105367.9611,148672.9455,
135343.5096,73087.3671,75604.4967,69973.8132,66367.3878,72191.2371,98143.1619,88506.7773,59722.086,59705.1591,54111.3165,49504.2126,50762.7774,68940.2766,67763.3592,50401.3383,50509.8696,51749.5161,
56774.814,56875.3797,65898.4131,58036.3659,45687.6945,44862.2592,41548.5696,37896.342,39405.8232,54220.8435,50035.9164,36221.5746,35067.5583,33099.0594,41550.561,33312.1392,42200.7531,36982.2894,26394.0156,
25939.9764,24163.6476,23641.9008,25495.8942,35859.1398,31292.8596,20783.2461,21233.3025,20770.302,20744.4138,19822.3956,31262.9886,28102.6368,17650.7739,17666.7051,17097.1647,16079.5593,17644.7997,26329.2951,
22879.1946,14800.0848,15000.2205,15227.2401,19307.6187,20345.1381,27219.4509,22791.573,14828.9601,26459.7318,65060.0337)
seq <- 1:length(val)
tmp <- data.frame(seq, val)

#log-normal function
fun <- function(x, mu, sd){
    a <- 1/(x * sqrt(2*pi*sd))
    b <- (log(x)-mu)^2 * -1
    c <- 2 * sd^2
    return (a * exp(b/c))
}

#fitting
library("minpack.lm")
fit <- nlsLM(val ~ a * fun(seq, b, c) + d,
             start = list(a=1, b=0.01, c=100, d=1), 
             control = nls.lm.control(maxiter = 1000),
             trace = TRUE,
             data=tmp)

#plotting
r <- sort(unique(tmp$seq))
xlim <- c(1,max(r))
ylim <- c(0, max(tmp$val) + max(tmp$val) * 0.3) 
plot(tmp$seq, tmp$val, cex=0.1, pch=19, xlim=xlim, ylim=ylim)
lines(r, predict(fit, data.frame(seq=r)), col="red")

8 pareto #

pareto <- function(y, m, s) {
  return (s * (1 + y/(m * (s-1)))^(-s-1)/(m * (s-1)))
}



install.packages("rmutil")
library("rmutil")

dpareto(y, m, s, log=FALSE)
ppareto(q, m, s)
qpareto(p, m, s)
hpareto(y, m, s)
rpareto(n, m, s)

9 覃燕 3螳 #

library("minpack.lm")
fit <- nlsLM(u ~ a*seq^b + c*seq^d + e*seq^f, 
             start = list(a=40000000, b=-0.05, c=8015503, d=-0.05, e=1674444, f=-0.05), 
             control = nls.lm.control(maxiter=1000),
             data=df, trace=T)
summary(fit)

10 譯手鍵 願鍵 #

calcPeriod <- function(x, y){
    incr <- x[2] - x[1]
    tmp <- spectrum(y, plot=FALSE)
    p <- (1/tmp$freq*incr)[which.max(tmp$spec)] # max of spectrum
    p
}

x <- seq(0, 50, by = 0.05)
y <- sin(x)
p <- calcPeriod(x, y) # result = 2pi

13 beta #

https://www.r-bloggers.com/fitting-distribution-x-to-data-from-distribution-y/
set.seed(3)
x <- rgamma(1e5, 2, .2)
plot(density(x))
 
# normalize the gamma so it's between 0 & 1
# .0001 added because having exactly 1 causes fail
xt <- x / ( max( x ) + .0001 )
 
# fit a beta distribution to xt
library( MASS )
fit.beta <- fitdistr( xt, "beta", start = list( shape1=2, shape2=5 ) )
 
x.beta <- rbeta(1e5,fit.beta$estimate[[1]],fit.beta$estimate[[2]])
 
## plot the pdfs on top of each other
plot(density(xt))
lines(density(x.beta), col="red" )
 
## plot the qqplots
qqplot(xt, x.beta)

蠍 蠍郁鍵..
企: : るジ讓曙 襦螻豺 企Ν 譯殊語. 襦螻豺
EditText : Print : Mobile : FindPage : DeletePage : LikePages : Powered by MoniWiki : Last modified 2019-04-12 13:15:32

螳讌襦 襯 觸願鍵 覲企 觜襦 覲願鍵 覲狩讌襷 訖襴襯 企Μ蠍 磯 螳 .