Corrected the split of the drift.R file.

Added example data IM1xIPC81.




git-svn-id: http://kriemhild.uft.uni-bremen.de/svn/drfit@61 d1b72e20-2ee0-0310-a1c4-ad5adbbefcdc

3 files changed, 287 insertions, 288 deletions

diff --git a/R/drdata.R b/R/drdata.R
new file mode 100644
index 0000000..9da7d96
--- /dev/null
+++ b/R/drdata.R
@@ -0,0 +1,39 @@
+drdata <- function(substances, experimentator = "%", db = "cytotox",
+    celltype="IPC-81",enzymetype="AChE",
+    organism="Vibrio fischeri",endpoint="Luminescence",whereClause="1",
+    ok="'ok','no fit'")
+{
+    library(RODBC) 
+    channel <- odbcConnect(db,uid="cytotox",pwd="cytotox",case="tolower")
+    slist <- paste(substances,collapse="','")
+    if (db == "cytotox") {
+        responsetype <- "viability"
+        testtype <- "celltype"
+        type <- celltype
+    } else {
+        if (db == "enzymes") {
+            responsetype <- "activity"
+            testtype <- "enzyme"
+            type <- enzymetype
+        } else {
+            responsetype <- "response"
+            testtype <- "organism"
+            type <- organism
+        }
+    }
+        
+    query <- paste("SELECT conc,",responsetype,",unit,experimentator,substance,",testtype,
+        ",ok FROM ", db, " WHERE substance IN ('",
+        slist,"') AND experimentator LIKE '",
+        experimentator,"' AND ",testtype," LIKE '",
+        type,"' AND ",
+        whereClause," AND ok in (",
+        ok,")",sep="")
+    if (db == "ecotox") query <- paste(query," AND type LIKE '",endpoint,"'",sep="")
+    data <- sqlQuery(channel,query)
+    odbcClose(channel)
+    names(data)[[1]] <- "dose"
+    names(data)[[2]] <- "response"
+    data$substance <- factor(data$substance,levels=substances)
+    return(data)
+}
diff --git a/R/drfit.R b/R/drfit.R
index 9da7d96..ff03d90 100644
--- a/R/drfit.R
+++ b/R/drfit.R
@@ -1,39 +1,253 @@
-drdata <- function(substances, experimentator = "%", db = "cytotox",
-    celltype="IPC-81",enzymetype="AChE",
-    organism="Vibrio fischeri",endpoint="Luminescence",whereClause="1",
-    ok="'ok','no fit'")
+drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
+        probit = TRUE, logit = FALSE, weibull = FALSE,
+        linlogit = FALSE, linlogitWrong = NA, allWrong = NA,
+        s0 = 0.5, b0 = 2, f0 = 0)
 {
-    library(RODBC) 
-    channel <- odbcConnect(db,uid="cytotox",pwd="cytotox",case="tolower")
-    slist <- paste(substances,collapse="','")
-    if (db == "cytotox") {
-        responsetype <- "viability"
-        testtype <- "celltype"
-        type <- celltype
-    } else {
-        if (db == "enzymes") {
-            responsetype <- "activity"
-            testtype <- "enzyme"
-            type <- enzymetype
+    if(!is.null(data$ok)) data <- subset(data,ok!="no fit") # Don't use data with 
+                                                            # ok set to "no fit"
+    substances <- levels(data$substance)
+
+    ri <- rix <- 0                        # ri is the index over the result rows
+                                          # rix is used later to check if any
+                                          # model result was appended
+    rsubstance <- array()                 # the substance names in the results
+    rndl <- vector()                      # number of dose levels
+    rn <- vector()                        # mean number of replicates 
+                                          # in each dose level
+    runit <- vector()                     # vector of units for each result row
+    rlhd <- rlld <- vector()              # highest and lowest doses tested
+    mtype <- array()                      # the modeltypes
+    sigma <- array()                      # the standard deviation of the residuals
+    logED50 <- vector()
+    stderrlogED50 <- vector()
+    a <- b <- c <- vector()
+
+    splitted <- split(data,data$substance)
+    for (i in substances) {
+        tmp <- splitted[[i]]
+        fit <- FALSE
+        if (length(tmp) != 0) {
+            unit <- levels(as.factor(as.vector(tmp$unit)))
+            cat("\n",i,": Fitting data...\n",sep="")
         } else {
-            responsetype <- "response"
-            testtype <- "organism"
-            type <- organism
+            unit <- ""
+            cat("\n",i,": No data\n",sep="")
         }
+        if (length(unit) > 1) {
+            cat("More than one unit for substance ",i,", halting\n\n",sep="")
+            break
+        }
+        if (length(tmp$response) == 0) {
+            nodata = TRUE
+        } else {
+            nodata = FALSE
+        }
+        rix <- ri
+        if (nodata) {
+            n <- ndl <- 0
+        } else {
+            ndl <- length(levels(factor(tmp$dose)))
+            n <- round(length(tmp$response)/ndl)
+            if (is.na(startlogED50[i])){
+                w <- 1/abs(tmp$response - 0.3)
+                startlogED50[[i]] <- sum(w * log10(tmp$dose))/sum(w)
+            }
+            highestdose <- max(tmp$dose)
+            lowestdose <- min(tmp$dose)
+            lhd <- log10(highestdose)
+            lld <- log10(lowestdose)
+            responseathighestdose <- mean(subset(tmp,dose==highestdose)$response)
+            responseatlowestdose <- mean(subset(tmp,dose==lowestdose)$response)
+            if (responseathighestdose < 0.5) {
+                inactive <- FALSE
+                if (responseatlowestdose < 0.5) {
+                    active <- TRUE
+                } else {
+                    active <- FALSE
+                    if (linlogit && 
+                        length(subset(linlogitWrong,linlogitWrong == i))==0 &&
+                        length(subset(allWrong,allWrong == i))==0) {
+                        m <- try(nls(response ~ linlogitf(dose,1,f,logED50,b),
+                                data=tmp,
+                                start=list(f=f0,logED50=startlogED50[[i]],b=b0)))
+                        if (!inherits(m, "try-error")) {
+                            fit <- TRUE
+                            ri <- ri + 1
+                            s <- summary(m)
+                            sigma[[ri]] <- s$sigma
+                            rsubstance[[ri]] <- i
+                            rndl[[ri]] <- ndl
+                            rn[[ri]] <- n
+                            runit[[ri]] <- unit
+                            rlld[[ri]] <- log10(lowestdose)
+                            rlhd[[ri]] <- log10(highestdose)
+                            logED50[[ri]] <- coef(m)[["logED50"]]
+                            if (logED50[[ri]] > rlhd[[ri]]) {
+                                mtype[[ri]] <- "no fit"
+                                logED50[[ri]] <- NA
+                                stderrlogED50[[ri]] <- NA
+                                a[[ri]] <- NA
+                                b[[ri]] <- NA
+                                c[[ri]] <- NA
+                            } else {
+                                mtype[[ri]] <- "linlogit"
+                                stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
+                                a[[ri]] <- coef(m)[["logED50"]]
+                                b[[ri]] <- coef(m)[["b"]]
+                                c[[ri]] <- coef(m)[["f"]]
+                            }
+                        }
+                    }
+
+                    if (probit &&
+                        length(subset(allWrong,allWrong == i))==0) {
+                        m <- try(nls(response ~ pnorm(-log10(dose),-logED50,scale),
+                                    data=tmp,
+                                    start=list(logED50=startlogED50[[i]],scale=1)))
+                        if (chooseone==FALSE || fit==FALSE) {
+                            if (!inherits(m, "try-error")) {
+                                fit <- TRUE
+                                ri <- ri + 1
+                                s <- summary(m)
+                                sigma[[ri]] <- s$sigma
+                                rsubstance[[ri]] <- i
+                                rndl[[ri]] <- ndl
+                                rn[[ri]] <- n
+                                runit[[ri]] <- unit
+                                rlld[[ri]] <- log10(lowestdose)
+                                rlhd[[ri]] <- log10(highestdose)
+                                logED50[[ri]] <- coef(m)[["logED50"]]
+                                c[[ri]] <- NA
+                                if (logED50[[ri]] > rlhd[[ri]]) {
+                                    mtype[[ri]] <- "no fit"
+                                    logED50[[ri]] <- NA
+                                    stderrlogED50[[ri]] <- NA
+                                    a[[ri]] <- NA
+                                    b[[ri]] <- NA
+                                } else {
+                                    mtype[[ri]] <- "probit"
+                                    stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
+                                    a[[ri]] <- coef(m)[["logED50"]]
+                                    b[[ri]] <- coef(m)[["scale"]]
+                                }
+                            }
+                        }
+                    }
+
+                    if (logit &&
+                        length(subset(allWrong,allWrong == i))==0) {
+                        m <- try(nls(response ~ plogis(-log10(dose),-logED50,scale),
+                                data=tmp,
+                                start=list(logED50=startlogED50[[i]],scale=1)))
+                        if (chooseone==FALSE || fit==FALSE) {
+                            if (!inherits(m, "try-error")) {
+                                fit <- TRUE
+                                ri <- ri + 1
+                                s <- summary(m)
+                                sigma[[ri]] <- s$sigma
+                                rsubstance[[ri]] <- i
+                                rndl[[ri]] <- ndl
+                                rn[[ri]] <- n
+                                runit[[ri]] <- unit
+                                rlld[[ri]] <- log10(lowestdose)
+                                rlhd[[ri]] <- log10(highestdose)
+                                logED50[[ri]] <- a[[ri]] <- coef(m)[["logED50"]]
+                                b[[ri]] <- coef(m)[["scale"]]
+                                c[[ri]] <- NA
+                                if (logED50[[ri]] > rlhd[[ri]]) {
+                                    mtype[[ri]] <- "no fit"
+                                    logED50[[ri]] <- NA
+                                    stderrlogED50[[ri]] <- NA
+                                    a[[ri]] <- NA
+                                    b[[ri]] <- NA
+                                } else {
+                                    mtype[[ri]] <- "logit"
+                                    stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
+                                }
+                            }
+                        }
+                    }
+
+                    if (weibull &&
+                        length(subset(allWrong,allWrong == i))==0) {
+                        m <- try(nls(response ~ pweibull(-log10(dose)+location,shape),
+                                data=tmp,
+                                start=list(location=startlogED50[[i]],shape=s0)))
+                        if (chooseone==FALSE || fit==FALSE) {
+                            if (!inherits(m, "try-error")) {
+                                fit <- TRUE
+                                ri <- ri + 1
+                                s <- summary(m)
+                                sigma[[ri]] <- s$sigma
+                                rsubstance[[ri]] <- i
+                                rndl[[ri]] <- ndl
+                                rn[[ri]] <- n
+                                runit[[ri]] <- unit
+                                rlld[[ri]] <- log10(lowestdose)
+                                rlhd[[ri]] <- log10(highestdose)
+                                a[[ri]] <- coef(m)[["location"]]
+                                b[[ri]] <- coef(m)[["shape"]]
+                                sqrdev <- function(logdose) {
+                                    (0.5 - pweibull( - logdose + a[[ri]], b[[ri]]))^2 
+                                }
+                                logED50[[ri]] <- nlm(sqrdev,startlogED50[[i]])$estimate
+                                c[[ri]] <- NA
+                                if (logED50[[ri]] > rlhd[[ri]]) {
+                                    mtype[[ri]] <- "no fit"
+                                    logED50[[ri]] <- NA
+                                    stderrlogED50[[ri]] <- NA
+                                    a[[ri]] <- NA
+                                    b[[ri]] <- NA
+                                } else {
+                                    mtype[[ri]] <- "weibull"
+                                    stderrlogED50[[ri]] <- NA
+                                }
+                            }
+                        }
+                    }
+
+                }
+
+            } else {
+                inactive <- TRUE
+            }
+        }
+        if (ri == rix) {          # if no entry was appended for this substance
+            ri <- ri + 1
+            rsubstance[[ri]] <- i
+            rndl[[ri]] <- ndl
+            rn[[ri]] <- n
+            if (nodata) {
+                rlld[[ri]] <- rlhd[[i]] <- NA
+                mtype[[ri]] <- "no data"
+                runit[[ri]] <- NA
+            } else {
+                rlld[[ri]] <- log10(lowestdose)
+                rlhd[[i]] <- log10(highestdose)
+                runit[[ri]] <- unit
+                if (inactive) {
+                    mtype[[ri]] <- "inactive"
+                } else {
+                    if (active) {
+                        mtype[[ri]] <- "active"
+                    } else {
+                        mtype[[ri]] <- "no fit"
+                    }
+                }
+            }
+            sigma[[ri]] <- NA
+            logED50[[ri]] <- NA
+            stderrlogED50[[ri]] <- NA
+            a[[ri]] <- NA
+            b[[ri]] <- NA
+            c[[ri]] <- NA
+        }
+    }
+    results <- data.frame(rsubstance, rndl, rn, rlld, rlhd, mtype, logED50, stderrlogED50, runit, sigma, a, b)
+    names(results) <- c("Substance","ndl","n","lld","lhd","mtype","logED50","std","unit","sigma","a","b")
+    if (linlogit) {
+        results$c <- c
     }
-        
-    query <- paste("SELECT conc,",responsetype,",unit,experimentator,substance,",testtype,
-        ",ok FROM ", db, " WHERE substance IN ('",
-        slist,"') AND experimentator LIKE '",
-        experimentator,"' AND ",testtype," LIKE '",
-        type,"' AND ",
-        whereClause," AND ok in (",
-        ok,")",sep="")
-    if (db == "ecotox") query <- paste(query," AND type LIKE '",endpoint,"'",sep="")
-    data <- sqlQuery(channel,query)
-    odbcClose(channel)
-    names(data)[[1]] <- "dose"
-    names(data)[[2]] <- "response"
-    data$substance <- factor(data$substance,levels=substances)
-    return(data)
+    rownames(results) <- 1:ri
+    return(results)
 }
diff --git a/R/linlogitf.R b/R/linlogitf.R
index 2af8b7e..2d69256 100644
--- a/R/linlogitf.R
+++ b/R/linlogitf.R
@@ -2,257 +2,3 @@ linlogitf <- function(x,k,f,mu,b)
 {
     k*(1 + f*x) / (1 + ((2*f*(10^mu) + 1) * ((x/(10^mu))^b)))
 }
-
-drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
-        probit = TRUE, logit = FALSE, weibull = FALSE,
-        linlogit = FALSE, linlogitWrong = NA, allWrong = NA,
-        s0 = 0.5, b0 = 2, f0 = 0)
-{
-    if(!is.null(data$ok)) data <- subset(data,ok!="no fit") # Don't use data with 
-                                                            # ok set to "no fit"
-    substances <- levels(data$substance)
-
-    ri <- rix <- 0                        # ri is the index over the result rows
-                                          # rix is used later to check if any
-                                          # model result was appended
-    rsubstance <- array()                 # the substance names in the results
-    rndl <- vector()                      # number of dose levels
-    rn <- vector()                        # mean number of replicates 
-                                          # in each dose level
-    runit <- vector()                     # vector of units for each result row
-    rlhd <- rlld <- vector()              # highest and lowest doses tested
-    mtype <- array()                      # the modeltypes
-    sigma <- array()                      # the standard deviation of the residuals
-    logED50 <- vector()
-    stderrlogED50 <- vector()
-    a <- b <- c <- vector()
-
-    splitted <- split(data,data$substance)
-    for (i in substances) {
-        tmp <- splitted[[i]]
-        fit <- FALSE
-        if (length(tmp) != 0) {
-            unit <- levels(as.factor(as.vector(tmp$unit)))
-            cat("\n",i,": Fitting data...\n",sep="")
-        } else {
-            unit <- ""
-            cat("\n",i,": No data\n",sep="")
-        }
-        if (length(unit) > 1) {
-            cat("More than one unit for substance ",i,", halting\n\n",sep="")
-            break
-        }
-        if (length(tmp$response) == 0) {
-            nodata = TRUE
-        } else {
-            nodata = FALSE
-        }
-        rix <- ri
-        if (nodata) {
-            n <- ndl <- 0
-        } else {
-            ndl <- length(levels(factor(tmp$dose)))
-            n <- round(length(tmp$response)/ndl)
-            if (is.na(startlogED50[i])){
-                w <- 1/abs(tmp$response - 0.3)
-                startlogED50[[i]] <- sum(w * log10(tmp$dose))/sum(w)
-            }
-            highestdose <- max(tmp$dose)
-            lowestdose <- min(tmp$dose)
-            lhd <- log10(highestdose)
-            lld <- log10(lowestdose)
-            responseathighestdose <- mean(subset(tmp,dose==highestdose)$response)
-            responseatlowestdose <- mean(subset(tmp,dose==lowestdose)$response)
-            if (responseathighestdose < 0.5) {
-                inactive <- FALSE
-                if (responseatlowestdose < 0.5) {
-                    active <- TRUE
-                } else {
-                    active <- FALSE
-                    if (linlogit && 
-                        length(subset(linlogitWrong,linlogitWrong == i))==0 &&
-                        length(subset(allWrong,allWrong == i))==0) {
-                        m <- try(nls(response ~ linlogitf(dose,1,f,logED50,b),
-                                data=tmp,
-                                start=list(f=f0,logED50=startlogED50[[i]],b=b0)))
-                        if (!inherits(m, "try-error")) {
-                            fit <- TRUE
-                            ri <- ri + 1
-                            s <- summary(m)
-                            sigma[[ri]] <- s$sigma
-                            rsubstance[[ri]] <- i
-                            rndl[[ri]] <- ndl
-                            rn[[ri]] <- n
-                            runit[[ri]] <- unit
-                            rlld[[ri]] <- log10(lowestdose)
-                            rlhd[[ri]] <- log10(highestdose)
-                            logED50[[ri]] <- coef(m)[["logED50"]]
-                            if (logED50[[ri]] > rlhd[[ri]]) {
-                                mtype[[ri]] <- "no fit"
-                                logED50[[ri]] <- NA
-                                stderrlogED50[[ri]] <- NA
-                                a[[ri]] <- NA
-                                b[[ri]] <- NA
-                                c[[ri]] <- NA
-                            } else {
-                                mtype[[ri]] <- "linlogit"
-                                stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
-                                a[[ri]] <- coef(m)[["logED50"]]
-                                b[[ri]] <- coef(m)[["b"]]
-                                c[[ri]] <- coef(m)[["f"]]
-                            }
-                        }
-                    }
-
-                    if (probit &&
-                        length(subset(allWrong,allWrong == i))==0) {
-                        m <- try(nls(response ~ pnorm(-log10(dose),-logED50,scale),
-                                    data=tmp,
-                                    start=list(logED50=startlogED50[[i]],scale=1)))
-                        if (chooseone==FALSE || fit==FALSE) {
-                            if (!inherits(m, "try-error")) {
-                                fit <- TRUE
-                                ri <- ri + 1
-                                s <- summary(m)
-                                sigma[[ri]] <- s$sigma
-                                rsubstance[[ri]] <- i
-                                rndl[[ri]] <- ndl
-                                rn[[ri]] <- n
-                                runit[[ri]] <- unit
-                                rlld[[ri]] <- log10(lowestdose)
-                                rlhd[[ri]] <- log10(highestdose)
-                                logED50[[ri]] <- coef(m)[["logED50"]]
-                                c[[ri]] <- NA
-                                if (logED50[[ri]] > rlhd[[ri]]) {
-                                    mtype[[ri]] <- "no fit"
-                                    logED50[[ri]] <- NA
-                                    stderrlogED50[[ri]] <- NA
-                                    a[[ri]] <- NA
-                                    b[[ri]] <- NA
-                                } else {
-                                    mtype[[ri]] <- "probit"
-                                    stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
-                                    a[[ri]] <- coef(m)[["logED50"]]
-                                    b[[ri]] <- coef(m)[["scale"]]
-                                }
-                            }
-                        }
-                    }
-
-                    if (logit &&
-                        length(subset(allWrong,allWrong == i))==0) {
-                        m <- try(nls(response ~ plogis(-log10(dose),-logED50,scale),
-                                data=tmp,
-                                start=list(logED50=startlogED50[[i]],scale=1)))
-                        if (chooseone==FALSE || fit==FALSE) {
-                            if (!inherits(m, "try-error")) {
-                                fit <- TRUE
-                                ri <- ri + 1
-                                s <- summary(m)
-                                sigma[[ri]] <- s$sigma
-                                rsubstance[[ri]] <- i
-                                rndl[[ri]] <- ndl
-                                rn[[ri]] <- n
-                                runit[[ri]] <- unit
-                                rlld[[ri]] <- log10(lowestdose)
-                                rlhd[[ri]] <- log10(highestdose)
-                                logED50[[ri]] <- a[[ri]] <- coef(m)[["logED50"]]
-                                b[[ri]] <- coef(m)[["scale"]]
-                                c[[ri]] <- NA
-                                if (logED50[[ri]] > rlhd[[ri]]) {
-                                    mtype[[ri]] <- "no fit"
-                                    logED50[[ri]] <- NA
-                                    stderrlogED50[[ri]] <- NA
-                                    a[[ri]] <- NA
-                                    b[[ri]] <- NA
-                                } else {
-                                    mtype[[ri]] <- "logit"
-                                    stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
-                                }
-                            }
-                        }
-                    }
-
-                    if (weibull &&
-                        length(subset(allWrong,allWrong == i))==0) {
-                        m <- try(nls(response ~ pweibull(-log10(dose)+location,shape),
-                                data=tmp,
-                                start=list(location=startlogED50[[i]],shape=s0)))
-                        if (chooseone==FALSE || fit==FALSE) {
-                            if (!inherits(m, "try-error")) {
-                                fit <- TRUE
-                                ri <- ri + 1
-                                s <- summary(m)
-                                sigma[[ri]] <- s$sigma
-                                rsubstance[[ri]] <- i
-                                rndl[[ri]] <- ndl
-                                rn[[ri]] <- n
-                                runit[[ri]] <- unit
-                                rlld[[ri]] <- log10(lowestdose)
-                                rlhd[[ri]] <- log10(highestdose)
-                                a[[ri]] <- coef(m)[["location"]]
-                                b[[ri]] <- coef(m)[["shape"]]
-                                sqrdev <- function(logdose) {
-                                    (0.5 - pweibull( - logdose + a[[ri]], b[[ri]]))^2 
-                                }
-                                logED50[[ri]] <- nlm(sqrdev,startlogED50[[i]])$estimate
-                                c[[ri]] <- NA
-                                if (logED50[[ri]] > rlhd[[ri]]) {
-                                    mtype[[ri]] <- "no fit"
-                                    logED50[[ri]] <- NA
-                                    stderrlogED50[[ri]] <- NA
-                                    a[[ri]] <- NA
-                                    b[[ri]] <- NA
-                                } else {
-                                    mtype[[ri]] <- "weibull"
-                                    stderrlogED50[[ri]] <- NA
-                                }
-                            }
-                        }
-                    }
-
-                }
-
-            } else {
-                inactive <- TRUE
-            }
-        }
-        if (ri == rix) {          # if no entry was appended for this substance
-            ri <- ri + 1
-            rsubstance[[ri]] <- i
-            rndl[[ri]] <- ndl
-            rn[[ri]] <- n
-            if (nodata) {
-                rlld[[ri]] <- rlhd[[i]] <- NA
-                mtype[[ri]] <- "no data"
-                runit[[ri]] <- NA
-            } else {
-                rlld[[ri]] <- log10(lowestdose)
-                rlhd[[i]] <- log10(highestdose)
-                runit[[ri]] <- unit
-                if (inactive) {
-                    mtype[[ri]] <- "inactive"
-                } else {
-                    if (active) {
-                        mtype[[ri]] <- "active"
-                    } else {
-                        mtype[[ri]] <- "no fit"
-                    }
-                }
-            }
-            sigma[[ri]] <- NA
-            logED50[[ri]] <- NA
-            stderrlogED50[[ri]] <- NA
-            a[[ri]] <- NA
-            b[[ri]] <- NA
-            c[[ri]] <- NA
-        }
-    }
-    results <- data.frame(rsubstance, rndl, rn, rlld, rlhd, mtype, logED50, stderrlogED50, runit, sigma, a, b)
-    names(results) <- c("Substance","ndl","n","lld","lhd","mtype","logED50","std","unit","sigma","a","b")
-    if (linlogit) {
-        results$c <- c
-    }
-    rownames(results) <- 1:ri
-    return(results)
-}