小编Luc*_*ini的帖子

我发现当我尝试使用绘图(...,add = T)覆盖多个栅格时,如果我尝试将多于3个栅格叠加在一起,则后续绘图不会正确对齐栅格.

我最初的意图是创建一个模拟的土地覆盖的分类地图,其中代表封面类的颜色的黑暗与我们的模型投影中的确定性不同.为此,我创建了一个简单的脚本,它将遍历每个封面类并使用从灰色(低确定性森林预测)到全覆盖颜色(例如,深绿色)的颜色渐变来绘制它(例如,森林,地图上的绿色)对于强烈预测的区域).我发现使用这种方法,在第3个封面添加到绘图之后,所有后续覆盖在绘图上的栅格都是任意错位的.我已经颠倒了封面类的绘图顺序,并且表现出相同的行为,这意味着它不是个别封面类栅格的问题.在Rstudio中更令人费解的是,当我使用变焦按钮仔细检查最终情节时,错位会恶化.

你对这种行为存在的原因有什么看法吗？最重要的是,您有任何建议的解决方案或解决方法吗？

下面链接中的代码和数据具有捕获的所有描述的行为. https://dl.dropboxusercontent.com/u/332961/r%20plot%20raster%20add%20issue.zip 转动plot_gradient = F看看你如何只是简单地将同一个栅格子集化并按顺序将子集添加到同一个图中可以复制这个问题.我已经尝试设置绘图设备图的范围(...,ext),但这不起作用.我也检查了每个封面光栅的范围是一样的.

下面是未对齐的封面类的图.绘制到jpeg设备将导致类似的图像(即,这不是Rstudio渲染的问题). 奇怪的是,如果我使用Rstudio放大图像,则不对齐是不同的 相比之下,这就是封面应该如何在景观中正确对齐

library(raster)
library(colorRamps)
raster_of_classes=raster("C:/r plot raster add issue/raster_of_classes.tif")
raster_of_certainty_of_classes=raster("C:/r plot raster add issue/raster_of_certainty_of_classes.tif")
endCols=c("darkorchid4", "darkorange3", "red3", "green4", "dodgerblue4") #colors to be used in gradients for each class
classes=unique(raster_of_classes)
minVal=cellStats(raster_of_certainty_of_classes, min)
tmp_i=1
addPlot=F
plot_gradient=F #this is for debug only
#classes=rev(classes) #turn this off and on to see how last 2 classes are mis aligned, regardless of plotting order
for (class in classes){
  raster_class=raster_of_classes==class #create mask for individual class
  raster_class[raster_class==0]=NA …

R中数据操作的一个常见任务是通过删除符合特定条件的行来对数据帧进行子类化.然而,在R中执行此操作的简单方法似乎在逻辑上不一致,甚至对无经验的人(如我自己)也很危险.

假设我们有一个数据框,我们想要排除属于"G1"处理的行:

Treatment=c("G1","G1","G1","G1","G1","G1","G2","G2","G2","G2","G2",
"G2","G3","G3","G3","G3","G3","G3")
Vals=c(runif(6),runif(6)+0.9,runif(6)-0.3)
data=data.frame(Treatment)
data=cbind(data, Vals)  

正如所料,下面的代码删除了与第一行标准匹配的数据帧行

to_del=which(data$Treatment=="G1")
new_data=data[-to_del,]
new_data

但是,与预期相反,使用此方法如果'which'命令找不到任何匹配行,则此代码将删除所有行,而不是将它们全部留下

to_del=which(data$Treatment=="G4")
new_data=data[-to_del,]
new_data

上面的代码导致数据框没有剩下的行,这没有任何意义(即,由于R没有找到符合我的删除标准的行,它删除了所有行).我的解决方法完成了这项工作,但我想如果没有所有这些条件语句,有一种更简单的方法可以做到这一点

###WORKAROUND
to_del=which(data$Treatment=="G4") #no G4 treatment in this particular data frame
if (length(to_del)>0){
  new_data=data[-to_del,]  
}else{
  new_data=data
}
new_data

有没有人有一个简单的方法来做到这一点即使没有符合指定条件的行？

很长一段时间以来,我一直在使用sfLapply来处理很多并行r脚本.然而,最近我已经深入研究并行计算,我一直在使用sfClusterApplyLB,如果单个实例不需要花费相同的时间来运行,那么可以节省大量时间.如果sfLapply将在加载新批处理之前等待批处理的每个实例完成(这可能导致空闲实例),完成任务的sfClusterApplyLB实例将立即分配给列表中的其余元素,因此可能会节省相当多的时间当实例没有花费相同的时间时.这让我质疑为什么我们在使用降雪时不想平衡我们的跑步？到目前为止我唯一发现的是,当并行脚本出现错误时,sfClusterApplyLB仍会在发出错误之前循环遍历整个列表,而sfLapply将在尝试第一批后停止.我还缺少什么？是否存在负载平衡的任何其他成本/缺点？下面是一个示例代码,显示了两者之间的差异

rm(list = ls()) #remove all past worksheet variables
working_dir="D:/temp/"
setwd(working_dir)
n_spp=16
spp_nmS=paste0("sp_",c(1:n_spp))
spp_nm=spp_nmS[1]
sp_parallel_run=function(sp_nm){
  sink(file(paste0(working_dir,sp_nm,"_log.txt"), open="wt"))#######NEW
  cat('\n', 'Started on ', date(), '\n') 
  ptm0 <- proc.time()
  jnk=round(runif(1)*8000000) #this is just a redundant script that takes an arbitrary amount of time to run
  jnk1=runif(jnk)
  for (i in 1:length(jnk1)){
    jnk1[i]=jnk[i]*runif(1)
  }
  ptm1=proc.time() - ptm0
  jnk=as.numeric(ptm1[3])
  cat('\n','It took ', jnk, "seconds to model", sp_nm)

  #stop sinks
  sink.reset <- function(){
    for(i in seq_len(sink.number())){
      sink(NULL)
    }
  }
  sink.reset()
}
require(snowfall)
cpucores=as.integer(Sys.getenv('NUMBER_OF_PROCESSORS'))

sfInit( parallel=T, cpus=cpucores) # …

好吧,这让我感到非常困惑和担心 - 作为例程的一部分,我一直在将变量的个别观察分类为TRUE或FALSE基于它们的值是高于还是低于/等于中值.但是,我在R中遇到的行为在执行这个简单的测试时很大程度上是出乎意料的.

所以采取这组观察:

data=c(0.6666667, 0.8333, 0.6666667, 0.8333, 0.8333, 0.75, 0.9999, 0.7499667, 0.25, 0.6666667, 0.1667, 0.7499667, 0.5, 0.2500333, 0.3333667, 0.0834, 0.0001, 0.2500333, 0.8333, 0.9999, 0.9999, 0.2500333, 0.2500333, 0.3333667, 0.9166, 0.5, 0.2500333, 0.4166667, 0.0001, 0.1667333, 0.6666333, 0.0834, 0.1667, 0.6666333, 0.9166, 0.1667, 0.7499333, 0.9166, 0.9166, 0.9166, 0.7499667, 0.7499667, 0.4166667, 0.5, 0.2500333, 0.9166, 0.6666667, 0.1667333, 0.25, 0.0001, 0.3333667, 0.0001, 0.25, 0.0834, 0.9999, 0.0834, 0.1667, 0.5, 0.2500333, 0.3333667, 0.9166, 0.9166, 0.8333, 0.9166, 0.75, 0.0834, 0.4166667, 0.5, 0.0001, 0.9999, 0.8333, 0.6666667, …

我有一个我正在创建的r代码,我想检测Windows中R的运行实例的数量,以便脚本可以选择是否运行特定的脚本集(即,如果已经有> 2个R的实例运行X,否则Y).

在R中有没有办法做到这一点？

编辑:这里有一些关于目的的信息:我有一套很长的脚本,用于在成千上万的情况下使用catnet库应用贝叶斯网络模型.此代码处理并输出每个案例的csv文件中的结果.我尝试过的大多数并行计算方案并不理想,因为它们抑制了许多内置的进度通知,因此我在R的不同实例上运行了一部分案例.我知道这有些过时,但它适用于我,所以我想要一种方法,根据运行的实例数自动使代码子集的数量.我现在手动通过在CMD中打开Rscript的多个实例来打开稍微不同配置的r文件来得到这样的东西:

cd "Y:\code\BN_code"
START "" "C:\Program Files\R\R-3.0.0\bin\x64\Rscript.exe" "process spp data3.r" /b
START "" "C:\Program Files\R\R-3.0.0\bin\x64\Rscript.exe" "process spp data3_T1.r" /b
START "" "C:\Program Files\R\R-3.0.0\bin\x64\Rscript.exe" "process spp data3_T2.r" /b
START "" "C:\Program Files\R\R-3.0.0\bin\x64\Rscript.exe" "process spp data3_T3.r" /b

EDIT2:

感谢下面的答案,这里是我在R中实现的'穷人的并行计算:所以如果你有任何长脚本必须应用于一长串案例列表,请使用下面的代码将长列表分成一个要提供给每个rscript实例的子列表数量:

#the cases that I need to apply my code to:
splist=c("sp01", "sp02", "sp03", "sp04", "sp05", "sp06", "sp07", "sp08", "sp09", "sp010", "sp11", "sp12", 
         "sp013", "sp014", "sp015", "sp16", "sp17", "sp018", "sp19", "sp20", "sp21", "sp22", "sp23", "sp24")
###automatic subsetting of cases based on …

在ggplot中,我们可以通过在数据框中指定具有条形高度的列来创建条形图

library("ggplot2")
library(plyr)
mm <- ddply(mtcars, "cyl", summarise, mmpg = mean(mpg))
ggplot(mm, aes(x = factor(cyl), y = mmpg)) + geom_bar(stat = "identity")

但是,我无法弄清楚如何制作一个类似的情节,指出了条形的顶部和底部.例如,使用下面的数据

df <- read.table(text = " id  min  max 
    Sp1     8.5          13.2     
 Sp2     11.7          14.5     
 Sp3     14.7          17.7     ", header=TRUE)

我们会得到一个非常类似的情节:

有什么建议？

在RI中,可以使用max/min命令轻松计算地理参考栅格堆栈中每个像元的最大/最小值.

set.seed(42)
require(raster)
r1 <- raster(nrows=10, ncols=10)
r2=r3=r4=r1
r1[]= runif(ncell(r1))
r2[]= runif(ncell(r1))+0.2
r3[]= runif(ncell(r1))-0.2
r4[]= runif(ncell(r1))
rs=stack(r1,r2,r3,r4)
plot(rs)
max(rs)
min(rs)

但是,我一直试图找到一种方法来找到堆栈中的第二高值.在我的例子中,堆栈上的每个栅格表示跨空间的特定模型的性能.我想比较第一个和第二个最佳值,以确定它的亚军最好的模型有多好,而不必将我的堆栈转换为矩阵然后再转换为栅格.任何想法或建议??

有没有一种有效的方法来计算重复指定值的向量部分的长度？例如,我想计算沿着每日降雨量值向量的无雨期的长度:

daily_rainfall=c(15, 2, 0, 0, 0, 3, 3, 0, 0, 10)

除了使用循环遍历向量的明显但笨重的方法之外,我可以通过哪种更清晰的方式获得所需的答案

rainless_period_length=c(3, 2)

鉴于上面的矢量？