小编gen*_*ser的帖子

我可以访问一个强大的大型集群.我是一个不错的R程序员,但对于shell命令(以及终端命令,除了使用ubuntu需要做的基本事情之外)都是全新的.

我想使用这个集群在R中运行一堆并行进程,然后我想将它们组合起来.具体来说,我遇到的问题类似于:

my.function <-function(data,otherdata,N){
    mod = lm(y~x, data=data)
    a = predict(mod,newdata = otherdata,se.fit=TRUE)
    b = rnorm(N,a$fit,a$se.fit)
    b
    }

r1 = my.function
r2 = my.function
r3 = my.function
r4 = my.function
...
r1000 = my.function

results = list(r1,r2,r3,r4, ... r1000)

以上只是一个愚蠢的例子,但基本上我想要并行执行1000次,然后对1000个进程的所有结果做一些事情.

如何同时向群集提交1000个作业,然后合并所有结果,如代码的最后一行？

我也欢迎任何有关编写RTFM的精心编写的手册/参考资料的建议.不幸的是,我发现的文件并不是特别容易理解.

提前致谢!

一个例子:

library(mgcv)
N=1000
x1 = seq(1:N)
x2 = log(x1)
x3 = sqrt(x1)
fac1 = ceiling(rnorm(N)*3)
fac2 = ceiling(runif(N)*3)
y = fac1*x2 + x1*x2 + x2 + x3*x2 + x2*(x1/x3)^(.8+fac2/10) + rnorm(N)*x2

mod = gam(y~
        s(as.factor(fac1),bs="re",by=x2) 
    +       s(x2)
    +   s(x1,by=x2)
    +   s(x3,by=x2)
    +   te(x1,x3,by=x2, by=as.factor(fac2))
    )   

最后一个张量不想让我两次互动.第一个是连续变量,它将模型矩阵中的每个项相乘,第二个是一个因子 - 它为每个因子级别创建不同的表面.

编辑:最后一个术语相当于te(x1*x2,x1*x3,by=as.factor(fac2)).但如果我发明了一个新变量x1x2 = x1*x2,我就失去了调用的能力predict.gam

我该如何编程呢？我需要上诉SmoothCon还是其他一些？如果是这样,一个如何实施的例子将非常有用.

谢谢!

(PS:我知道模型的异方性.特征,而不是bug.)

我正在将随机森林训练为多级数据，本质上将其视为非参数回归模型。我将未观察到的组级异质性视为随机森林训练过程之外的修正。但在随机森林训练过程中，我希望每棵树都是从横截面单元的随机引导样本（或子样本）而不是观察中生长出来的。因此，假设我的数据是对许多个体的多次观察，我想引导个体，而不是对个体的观察。

下面的虚拟示例表明这并strata不能解决我的问题。

> N <- 1000
> p <- 100
> A <- matrix(rnorm(p^2),p)
> library(MASS)
> X <- mvrnorm(N, rep(0,p), A %*%t(A))
> B <- rnorm(p)
> fac <- sample(1:1000 %% 10 +1)
> y <- log(fac + exp(X%*%B)^{1/fac}) + rnorm(N, sd = 10)
> fac <- as.factor(fac)
> library(randomForest)
> forest <- randomForest(y = y, x = cbind(X, fac), ntree = 1, keep.inbag = TRUE, replace = FALSE
+  , strata = fac                        #Stratify …

我试图让箱形图从具有最低平均值的因子发展到具有最高平均值的因子.这是一个简单的例子:

a = rnorm(10,mean=3,sd=4)
b = rnorm(10,mean=-1,sd=2)
c = rnorm(10,mean=5,sd=6)
d = rnorm(10,mean=-3,sd=1)
e = rnorm(10,mean=0,sd=.5)

labs = c(rep("a",10),rep("b",10),rep("c",10),rep("d",10),rep("e",10))
mean =     c(rep(mean(a),10),rep(mean(b),10),rep(mean(c),10),rep(mean(d),10),rep(mean(e),10))
data = c(a,b,c,d,e)
df = data.frame(labs,data,mean)
df = df[order(df$mean),]
boxplot(data~labs,data=df)
#They are not ordered
df$labs = ordered(df$labs, levels=levels(df$labs))
boxplot(data~labs,data=df)
#It doesn't work

我怎样才能得到左边最小的因子,随着我向右走？这有几个主题,但他们的方法对我不起作用.(也许是因为我的数据格式？)

奖励积分帮助我将x轴上的字母旋转180度.

提前致谢!

我正在尝试从此处获取全球净初级生产力的时间序列估计值，将其下载到基于美国的远程计算机上（我在肯尼亚，手机上网速度很慢），并将其从当前格式转换为ascii，使用如下函数粗化它们：

coarsen = function(input,size_out){
    if (dim(input)[1] != dim(input)[2]) {print("Matrix not square.")} else {    
    dimsize = dim(input)[1]
    fac = dimsize/size_out
    x = suppressWarnings(matrix(c(rep(1, fac), rep(0, dimsize)), dimsize, dimsize/fac))
    t(x) %*% input %*% x}
    }

然后下载它们，以便我可以与他们合作。使用R而不是专用的GIS程序可以做到吗？如果是这样，我该怎么办？首选R解决方案-我将使用无法在其上安装程序的群集，但可以安装R软件包。

我主要研究观察数据，但我阅读了很多实验性硬科学论文，这些论文以方差分析表的形式报告结果，用字母表示组之间差异的显着性，然后是 f 的 p 值-stat 表示本质上是因子变量回归的联合显着性。这是我从谷歌图片搜索中提取的一个例子。

我认为，在我继续尝试以各种方式控制它们之前，这可能是在观察数据集中呈现有关分组差异（或缺乏分组）的汇总统计数据的有用方法。我不确定这些字母通常代表什么测试（Tukey 什么？），但成对 t 检验很适合我的目的。

我的主要问题是：如何从 R 中的因子变量回归中获得这样的输出，以及如何将其无缝导出到乳胶中？

以下是一些示例数据：

var = c(3500,600,12400,6400,1500,0,4400,400,900,2000,350,0,5800,0,12800,1200,350,800,2500,2000,0,3200,1100,0,0,0,0,0,1000,0,0,0,0,0,12400,6000,1700,3500,3000,1000,0,0,3500,5000,1000,3600,1600,3500,0,900,4200,0,0,0,0,1700,1250,500,950,500,600,1080,500,980,670,1200,600,550,4000,600,2800,650,0,3700,12500,0,0,0,1200,2700,0,NA,0,0,0,3700,2000,3500,0,0,0,3500,800,1400,0,500,7000,3500,0,0,0,0,2700,0,0,0,0,2000,5000,0,0,7000,0,4800,0,0,0,0,1800,0,2500,1600,4600,0,2000,5400,4500,3200,0,12200,0,3500,0,0,2800,3600,3000,0,3150,0,0,3750,2800,0,1000,1500,6000,3090,2800,600,0,0,1000,3800,3000,0,800,600,1200,0,240,1000,300,3600,0,1200,300,2700,NA,1300,1200,1400,4600,3200,750,300,750,1200,700,870,900,3200,1300,1500,1200,0,960,1800,8000,1200,NA,0,1080,1300,1080,900,700,5000,1500,3750,0,1400,900,1400,400,3900,0,1400,1600,960,1200,2600,420,3400,2500,500,4000,0,4250,570,600,4550,2000,0,0,4300,2000,0,0,0,0,NA,0,2060,2600,1600,1800,3000,900,0,0,3200,0,1500,3000,0,3700,6000,0,0,1250,1200,12800,0,1000,1100,0,950,2500,800,3000,3600,3600,1500,0,0,3600,800,0,1000,1600,1700,0,3500,3700,3000,350,700,3500,0,0,0,0,1500,0,400,0,0,0,0,0,0,0,500,0,0,0,0,5600,0,0,0)
factor = as.factor(c(5,2,5,5,5,3,4,5,5,5,3,1,1,1,5,3,6,6,6,5,5,5,3,5,3,3,3,3,4,3,3,3,4,3,5,5,3,5,3,3,3,3,5,3,3,3,3,3,5,5,5,5,5,3,3,5,3,5,5,3,5,5,4,3,5,5,5,5,5,5,4,5,3,5,4,4,3,4,3,5,3,3,5,5,5,3,5,5,4,3,3,5,5,4,3,3,5,3,3,4,3,3,3,3,5,5,3,5,5,3,3,5,4,3,3,3,4,4,5,3,1,5,5,1,5,5,5,3,3,4,5,5,5,3,3,4,5,4,5,3,5,5,5,3,3,3,3,3,3,3,3,3,3,3,4,3,3,3,3,3,3,3,4,5,4,6,4,3,5,5,3,5,3,3,4,3,5,5,5,3,5,3,3,5,5,5,3,4,3,3,3,5,3,5,3,5,5,3,5,3,5,5,5,5,5,3,5,3,5,3,4,5,5,5,6,5,5,5,5,4,5,3,5,3,3,5,4,3,5,3,4,5,3,5,3,5,3,1,5,1,5,3,5,5,5,3,6,3,5,3,5,2,5,5,5,1,5,5,6,5,4,5,4,3,3,3,5,3,3,3,3,5,3,3,3,3,3,3,5,5,5,4,4,4,5,5,3,5,4,5,5,4,3,3,3,4,3,5,5,4,3,3))

对它们做一个简单的回归，你会得到以下结果

m = lm((var-mean(var,na.rm=TRUE))~factor-1)
summary(m)
Call:
lm(formula = (var - mean(var, na.rm = TRUE)) ~ factor - 1)

Residuals:
    Min      1Q  Median      3Q     Max 
-2040.5 -1240.2  -765.5   957.1 10932.8 

Coefficients:
        Estimate Std. Error t value Pr(>|t|)  
factor1   -82.42     800.42  -0.103   0.9181  
factor2  -732.42    1600.84  -0.458   0.6476  
factor3  -392.17     204.97  -1.913   0.0567 .
factor4   -65.19     377.32  -0.173   0.8629  
factor5 …

我正在编写一个函数,我想为用户提供并行运行的选项.foreach出于各种原因,我更喜欢竞争对手.到目前为止,我已经将其设置为类似于下面的虚拟函数:

library(foreach)
myfun <- function(parallel = TRUE){
  if (parallel == TRUE){
    require(doMC)
    registerDoMC(detectCores())
    foreach(i = 1:10) %dopar% {
      print(i)
    }
  } else {
    foreach(i = 1:10) %do% {
      print(i)
    }
  }
  return('OK!')
}
myfun()

这比它看起来需要的时间更长,并且冒着我可能不记得反映对底部顶部所做的任何更改的风险.我宁愿做类似以下的事情(非工作),但我不知道是否可以这样做:

myfun <- function(parallel = TRUE){
  if (parallel == TRUE){
    require(doMC)
    registerDoMC(detectCores())
  }
  foreach(i = 1:10) ifelse(parallel == TRUE, %dopar%, %do%) {
    print(i)
  }
  return('OK!')
}

任何选择%do%vs的方法都%dopar%取决于函数的arg？我想我可以只注册1个内核并且%dopar%每次调用,但是这会引发一个关于没有注册并行后端的警告,这可能会让那些不是我的人感到困惑.

"y"和"g"的底部部分在此地图中被截断.在大豆下看"旱地"和"灌溉".

什么设置控制这个？这是我的theme_map功能

theme_map <- function (base_size = 12, base_family = "") {
  theme_gray(base_size = base_size, base_family = base_family) %+replace% 
    theme(
      axis.line=element_blank(),
      axis.text.x=element_blank(),
      axis.text.y=element_blank(),
      axis.ticks=element_blank(),
      axis.ticks.length=unit(0.3, "lines"),
      axis.ticks.margin=unit(0.5, "lines"),
      axis.title.x=element_blank(),
      axis.title.y=element_blank(),
      legend.background=element_rect(fill="white", colour=NA),
      legend.key=element_rect(colour="white"),
      legend.key.size=unit(1.2, "lines"),
      legend.position="right",
      legend.text=element_text(size=rel(0.8)),
      legend.title=element_text(size=rel(0.8), face="bold", hjust=0),
      panel.background=element_blank(),
      panel.border=element_blank(),
      panel.grid.major=element_blank(),
      panel.grid.minor=element_blank(),
      panel.spacing=unit(0, "lines"),
      plot.background=element_blank(),
      plot.margin=unit(c(1, 1, 0.5, 0.5), "lines"),
      plot.title=element_text(size=rel(1.2)),
      strip.background=element_rect(fill="white", colour="white"),
      strip.text.x=element_text(size=rel(.8)),
      strip.text.y=element_text(size=rel(0.8), angle=-90) 
    )   
}

ggplot2在花了很长时间作为基础图形的忠实奉献者之后,我只是蘸了我的脚趾.我theme_map在互联网上找到了这个功能,所以我还不明白.有人可以告诉我,我需要调整哪些设置来解决我的问题？

这是我生成(其中一个)这些地图的代码,但我假设相关部分在theme_map函数中

mapcorn <- ggplot(USA_fort_premade) + 
    aes(long,lat,group=group) + 
    geom_polygon(aes(x=long,y=lat, group=group, fill=yield.mean), 
                 data=REAP_fort_corn)+
    geom_path(color="grey", …

我正在修改Brieman的随机森林程序(我不知道C/C++),所以我在R中从头开始编写自己的RF变体.我的程序和标准程序之间的区别基本上只是在如何计算终端节点中的分割点和值 - 一旦我在森林中有一棵树,它可以被认为与典型的RF算法中的树非常相似.

我的问题是它的预测速度很慢,而且我很难想办法让它更快.

测试树对象链接在这里,和一些测试数据被链接在这里.您可以直接下载,也可以在下面加载它(如果已repmis安装).他们被称为testtree和sampx.

library(repmis)
testtree <- source_DropboxData(file = "testtree", key = "sfbmojc394cnae8")
sampx <- source_DropboxData(file = "sampx", key = "r9imf317hpflpsx")

编辑:不知怎的,我还没有真正学习如何使用github.我上传所需的文件到存储库在这里 -我无法弄清楚如何在此刻得到一个永久的歉意......

它看起来像这样(使用我编写的绘图函数):

这里有一些关于对象的结构:

1> summary(testtree)
         Length Class      Mode   
nodes       7   -none-     list   
minsplit    1   -none-     numeric
X          29   data.frame list   
y        6719   -none-     numeric
weights  6719   -none-     numeric
oob      2158   -none-     numeric
1> summary(testtree$nodes)
     Length Class  Mode
[1,] 4      -none- list
[2,] 8      -none- list
[3,] …

我有如下所示的数据框：

import time
import pandas as pd
import numpy as np

N = 3
l = []
for i in range(N):
    n = np.random.choice(5)+2
    l += [pd.DataFrame(dict(ID = np.repeat(i, n),
                            t = list(range(n)),
                            X = np.random.normal(size = n)))]
df = pd.concat(l)

df
Out[85]: 
   ID  t         X
0   0  0  0.992300
1   0  1  0.226487
2   0  2 -0.731178
3   0  3  0.748376
4   0  4  1.269106
0   1  0  0.512957
1   1  1 -1.274963
2   1  2  0.186314
3   1 …