小编den*_*nis的帖子

I am facing a problem I do not manage to solve. I would like to use nlme or nlmODE to perform a non linear regression with random effect using as a model the solution of a second order differential equation with fixed coefficients (a damped oscillator).

I manage to use nlme with simple models, but it seems that the use of deSolve to generate the solution of the differential equation causes a problem. Below an example, and the problems I …

我试图使用nlme和拟合一阶微分模型lsoda。这是基本思想：我首先定义允许生成微分方程解的函数：

library(deSolve)

ODE1 <- function(time, x, parms) {with(as.list(c(parms, x)), {
  import <- excfunc(time)
  dS <- import*k/tau - (S-yo)/tau 
  res <- c(dS)
  list(res)})}


solution_ODE1 = function(tau1,k1,yo1,excitation,time){
  excfunc <- approxfun(time, excitation, rule = 2)
  parms  <- c(tau = tau1, k = k1, yo = yo1, excfunc = excfunc)
  xstart = c(S = yo1)
  out <-  lsoda(xstart, time, ODE1, parms)
  return(out[,2])
}

然后，根据两个ID的公式生成数据：

time <- 0:49
excitation <- c(rep(0,10),rep(1,10),rep(0,10),rep(1,10),rep(0,10))
simu_data <- data.frame(signal = c(solution_ODE1(3,2,0.1,excitation,time)+rnorm(length(time),0,0.1),
                                   solution_ODE1(3.2,1.5,0.3,excitation,time)+rnorm(length(time),0,0.1)),
                        time = rep(time,2),
                        excitation …

几周前我发现ggforce，它有一个很好的绘制椭圆的功能。但我无法在日志图中使用它。下面是一个例子：

我想用椭圆圈出这个组

library(ggforce)
library(ggplot2)

ggplot(mtcars)+
  geom_point(aes(hp,disp))+
  geom_ellipse(aes(x0 = 230, y0 = 450, a = 80, b = 30, angle = -10))

但我想在对数图中这样做。如果我天真地做

ggplot(mtcars)+
  geom_point(aes(hp,disp))+
  geom_ellipse(aes(x0 = 230, y0 = 450, a = 80, b = 30, angle = -10))+
  scale_y_log10()

我得到一个巨大的椭圆：

看起来椭圆参数没有对数转换。我可以尝试减少参数轴以获得对数轴上的合适大小，例如：

ggplot(mtcars)+
  geom_point(aes(hp,disp))+
  scale_y_log10()+
  geom_ellipse(aes(x0 = 230, y0 = 450, a = 80, b = 0.05, angle =0))

哪个有效：

但只有当角度为0时，如果不是，两个wxis是混合的，我无法得到我想要的椭圆：

ggplot(mtcars)+
  geom_point(aes(hp,disp))+
  scale_y_log10()+
  geom_ellipse(aes(x0 = 230, y0 = 450, a = 80, b = 0.05, angle = …

我想使用 h2o 进行Rglm 回归，但具有随机效应（HGLM，从本页看来是可能的）。我还没有设法让它工作，并出现我不明白的错误。

这是我的工作示例：我用辛普森悖论定义了一个数据集：全球呈上升趋势，但每组呈下降趋势

library(tidyverse)
library(ggplot2)
library(h2o)
library(data.table)

global_slope <- 1
global_int <- 1

Npoints_per_group <- 50
N_groups <- 10
pentes <- rnorm(N_groups,-1,.5)

centers_x <- seq(0,10,length = N_groups)
center_y <- global_slope*centers_x + global_int

group_spread <- 2

group_names <- sample(LETTERS,N_groups)

df <- lapply(1:N_groups,function(i){
  x <- seq(centers_x[i]-group_spread/2,centers_x[i]+group_spread/2,length = Npoints_per_group)
  y <- pentes[i]*(x- centers_x[i])+center_y[i]+rnorm(Npoints_per_group)
  data.table(x = x,y = y,ID = group_names[i])
}) %>% rbindlist()

您可以识别出类似于辛普森悖论 wiki 页面示例的内容：

ggplot(df,aes(x,y,color = as.factor(ID)))+
  geom_point()

没有随机效应的线性回归呈现增加趋势：

lm(y~x,data = df) %>% …

我有一个问题,我无法在data.table中正确解决.我有以下数据:

plouf <- data.table(  ID = rep(LETTERS[1:10],each = 10) )
plouf[,c(paste0("X",1:10)) := lapply(1:10,function(x){sample(10,100,replace = T)})]

有两件事阻碍了我的时间:

 col <- "X1"
 plouf[get(col) > 5, .(col = get(col)[1]) ,by = ID]
    ID col
 1:  A   7
 2:  B   7
 3:  C   9
 4:  D   6
 5:  E   8
 6:  F   7
 7:  G   6
 8:  H   7
 9:  I   6
10:  J   7

该列名为"col"而不是"X1".我试过了eval,get没有得到它.

同样的:

 col <- 1
 plouf[get(paste0("X",col)) > 5, .(paste0("X",col) = get(paste0("X",col))[1]) ,by = ID]

Error: unexpected …

当我coxph(Surv(Time, Status)~Class, data = df)在如下数据集上运行时，它总是按字母顺序设置参考组，在这种情况下，MutantA 将是参考组。有没有办法告诉它让 WT 成为参考组？

df

Time    Status  Class
3   1   WT
4   1   WT
5   1   WT
7   1   WT
7   1   WT
7   1   WT
7   1   WT
2   1   WT
2   1   WT
2   1   WT
5   1   WT
6   1   WT
7   1   WT
8   1   MutantA
9   1   MutantA
2   1   MutantA
12  1   MutantA
3   1   MutantA
4   1   MutantA
5   1   MutantA
7   1   MutantA
7   1   MutantA …

我有以下问题 - 在这样开始的脚本中：

modules::import("modules")
modules::import("futile.logger")
modules::import("data.table")
modules::import("REDCapR")

作为在其他脚本中使用的模块modules::use，我想使用redcap_write包函数的修改版本REDCapR。

我不知道如何继续。在我看来，有两种可能性：

• 使用本地存储的修改版本redcap_write。那太好了，因为这将是一个很容易分享的修改。但我不知道如何强制R用redcap_write我本地修改的版本替换包的功能。modules::use只会导入修改后的函数，但不会替换包版本redcap_write

• REDCapR安装我在此处创建的软件包的分叉版本https://github.com/dmongin/REDCapR/tree/overwrite。但我不知道如何以简单的方式进行（卸载REDCapR以安装我的分叉版本会有点麻烦：我们在各个开发人员/用户之间共享代码，每个人都应该卸载并重新安装软件包）

我在绘图方面遇到了困难：我想删除ggplot绘图中填充图例的一部分，同时保持自动着色。这是一个例子：

library(ggplot2)

df1 <- data.frame(x = 1:20,y1 = rnorm(20,2,0.2),y2 = sqrt(1:20))
df2 <- data.frame(x1 = c(1,5,10),x2 = c(5,10,20),color2 = as.factor(1:3))

ggplot(data=df1) +
  geom_rect(data = df2,
            aes(xmin = x1,
                xmax = x2,
                ymin = 0,
                ymax = Inf,
                fill = color2),
            color = "black",
            size = 0.3,
            alpha = 0.2)+
  geom_bar(aes(x = x,
               y= y1,
               fill = "daily"),
           stat='identity',
           width = 0.75,
           size = 0.1,
           alpha = 0.5) +
  geom_line(aes(x = x,
                y =y2,
                color = "somthing"),
            size = …

我想在区间[-50,50]中生成100个正态分布的随机数.但是在下面的代码中,生成的随机数范围是[-50,50].

n <- rnorm(100, -50,50)
plot(n)