在R中计算处理后的差异

Question

我对 R 中的面板数据有疑问。

我的数据基本上是这样的：

Year  Name       Variable    Treatment
2000  CompanyA   10          0
2001  CompanyA   10          0
2002  CompanyA   10          1
2003  CompanyA   10          0
2004  CompanyA   12          0
2005  CompanyA   12          0
1999  CompanyB    5          1
2000  CompanyB    5          1
2001  CompanyB    5          0
2002  CompanyB    5          0
2003  CompanyB    6          0
2004  CompanyB    5          0
2005  CompanyB    6          0
2006  CompanyB    6          0

有没有机会在R中计算治疗前后因变量的差异（关于一定的时滞）？

不幸的是，我只有不平衡的面板数据。计算的目的是从中得出一个虚拟变量。这将显示两年后因变量是否增长。然后，我想对其进行 clogit 回归。

编辑

我需要知道因变量在治疗后是否发生了变化。 所以我需要某种代码来为我的变量的每一个积极变化计算一个虚拟变量。

输出应该是这样的：

Year  Name       Variable    Treatment   Dummy
2000  CompanyA   10          0           0
2001  CompanyA   10          0           0
2002  CompanyA   10          1           0
2003  CompanyA   10          0           0
2004  CompanyA   12          0           1
2005  CompanyA   12          0           1
1999  CompanyB    5          1           0
2000  CompanyB    5          1           0
2001  CompanyB    5          0           0
2002  CompanyB    5          0           0
2003  CompanyB    6          0           1
2004  CompanyB    5          0           0
2005  CompanyB    6          0           0
2006  CompanyB    6          0           0

所以我可以对其进行条件 logit 回归，并将治疗（包括其他变量）与在一定时间滞后后对我的因变量的积极影响联系起来。

Answer 1

根据评论中的说明更新了答案；除了简单的比较（开/关处理，A 部分）之外，我还按照要求（B 部分）结合了时间进程的方法。
请注意，在许多情况下，代码需要适应确切的问题（如何处理成为治疗阴性，然后甚至可能再次成为 pos 的tose？从开始（或停止后）预测治疗效果的持续时间是多少？ o ftretment？这些问题可能更像是一个概念问题而不是 R 问题，但我试图提供一些如何实现这些问题的起点。

#### sample data (added and changed some data to demonstarte sorting of the years ####
# and pos Treatment at first time point):

text <- "Year  Name       Variable    Treatment
2000  CompanyA   10          0
2001  CompanyA   10          0
2002  CompanyA   10          1
2003  CompanyA   10          0
2004  CompanyA   12          0
2010  CompanyA   15          1
2005  CompanyA   12          0
1999  CompanyB    5          0
2000  CompanyB    5          1
2001  CompanyB    5          0
2002  CompanyB    5          0
2003  CompanyB    6          0
2004  CompanyB    5          0
2005  CompanyB    6          0
2006  CompanyB    6          0
2001  CompanyC    5          1
2006  CompanyC    9          1"

df <- read.table(text=text, header=TRUE)
str(df)
head(df)

#### A) Simple way: just compare on/off treatment subject ####

mean(df[df$Treatment==1, "Variable"]) - mean(df[df$Treatment==0, "Variable"]) 


#### B) Compare within each company, take into consideration also the time course ####

# split to list according to company names, to analyse them separately
Name.u <- as.character(unique(df$Name))  # unique Company names
L <- sapply(Name.u, function(n) df[df$Name==n, ], simplify=FALSE)             
str(L)
L  # a list of dataframes, one dataframe for each company

## deal with special cases that may influence the concept of theanalysis
# sort for year (assuming there are nor ties)
L <- sapply(Name.u, function(n) L[[n]][order(L[[n]]$Year), ], simplify=FALSE) 
# posibly ignore those who were already treatet at study entry already
L.del <- sapply(Name.u, function(n) ifelse(L[[n]][1, "Treatment"]==1, TRUE, FALSE), simplify=TRUE) 
L[L.del] <- NULL
Name.u <- Name.u[!L.del]
str(L); L # note that CompanyC was deleted because of Treatment==1 at start

## display treatment duration etc.
LL <- function(L.n) {
  L.n$diff <- c(0, diff(L.n$Treatment))
  # stopifnot(sum(L.n$diff!=0) == 1)   # more than one status change - need clarification how this should be handled, see also lines below
  # ALL status change to "treated" (possibly more than one!)
  Rx.start <- which(L.n$diff==1) 
  # duration since FIRST documented treatment
  L.n$RxDurSinceFirst <- L.n$Year - min(L.n$Year[Rx.start])  
  L.n$RxDurReal <- L.n$RxDur
  # need to define what to do with those who are Treatment negative at THIS  time ...
  L.n$RxDurReal[L.n$Treatment==0] <- NA   
  # ... and those who became Treatment neg before or now
  L.n$RxDurReal[sapply(1:nrow(L.n), function(row.i) row.i >= min(which(L.n$diff==-1)))] <- NA  
  return(L.n)
}
str(LL)

# L2 is a new list of the same structure as L, but with more information 
# (more columns in each dataframe element)
L2 <- sapply(Name.u, function(n) LL(L[[n]]), simplify=FALSE)
str(L2)
L2

# for a company n one can then do (and of course further vectorize):
n <- Name.u[1]
str(L2[[n]])
L2[[n]]

# for a company n one can then compare RxDurSinceFirst, RxDurReal or 
# whateveryou want (and of course further vectorize):
(Var.before <- L2[[n]]$Variable[ L2[[n]]$RxDurSinceFirst <  0 ] )
(Var.after  <- L2[[n]]$Variable[ L2[[n]]$RxDurSinceFirst >= 0 ] )
t.test(Var.before, Var.after)  # works of course only if enough observations

# or on/off Treatment within one group, and use the means of each group 
# for further paired t.test/ U-test etc.
(Var.OnRx  <- L2[[n]]$Variable[ L2[[n]]$Treatment ==  0 ] )
(Var.OffRx <- L2[[n]]$Variable[ L2[[n]]$Treatment ==  1 ] )

### End ###

Answer 2

这是一个我认为会让您非常接近的答案。我的代码突出显示了处理前变量的任何变化。请注意，这不是最优雅的代码，或多或少是草稿版本，但我必须打包，我认为这可能仍然有用。

首先，这是您的表格的 dput。只需运行它即可加载表格。

dfx <- structure(list(Year = c(2000L, 2001L, 2002L, 2003L, 2004L, 2005L, 
1999L, 2000L, 2001L, 2002L, 2003L, 2004L, 2005L, 2006L), Name = c("CompanyA", 
"CompanyA", "CompanyA", "CompanyA", "CompanyA", "CompanyA", "CompanyB", 
"CompanyB", "CompanyB", "CompanyB", "CompanyB", "CompanyB", "CompanyB", 
"CompanyB"), Variable = c(10L, 10L, 10L, 10L, 12L, 12L, 5L, 5L, 
5L, 5L, 6L, 5L, 6L, 6L), Treatment = c(0L, 0L, 1L, 0L, 0L, 0L, 
1L, 1L, 0L, 0L, 0L, 0L, 0L, 0L), Dummy = c(0L, 0L, 0L, 0L, 1L, 
1L, 0L, 0L, 0L, 0L, 1L, 0L, 0L, 0L)), .Names = c("Year", "Name", 
"Variable", "Treatment", "Dummy"), class = "data.frame", row.names = c(NA, 
-14L))

然后我创建了一个辅助变量 (has_treatment)，用于说明某年（行）是否接受过治疗。这是此函数中的前两行。

然后是一个简单的条件语句，我在其中测试一个案例是否已经过治疗，以及变量是否与治疗前的变量不同。

foo <- function(dfx){
      dfx[(Position( isTRUE, diff(dfx$Treatment) == -1)+1)  : nrow(dfx), "has_treatment" ] <- 1 

      dfx[1:(Position( isTRUE, diff(dfx$Treatment) == -1))  , "has_treatment" ] <- 0 

      dfx[dfx$has_treatment == 1 & 
              ((dfx[dfx$Treatment == 1, "Variable"] == 
                  dfx[, "Variable"])==FALSE) ,"dummy"] <- 1
  return(dfx)
}

然后我在 ddply 中运行它。如果您不熟悉 ddply 和 plyr 包，我强烈建议您学习一下。

library(plyr)

ddply(test, .variables = "Name", foo   )

同样，这并不是您想要的，但原则上它应该让您走上正轨。我会尝试再试一次，但我必须跑。

另外，正如一些人可能会评论的那样，这不是最优雅的方式，并且可能有更快、更有效的方式。

无论如何，我希望它会有所帮助。

Answer 3

或者，

diff(by(df$Variable, df$Treatment, FUN=mean))
#[1] -1.242424