Hoas

  • 3 天前
  • 注册于 2018年11月5日

    AgnesZhang 好久之前的了,我现在也好久没用了,忘了是哪些包。

  • 如题,偏态分布,想用《概率论及数理统计》书中的公式,但是不会实践,给个例子吧:

    library(lmom)
    
    x <- quape3(runif(100), para = c(7.5, 7.5 * 0.45, 0.45 * 3))
    z <- quape3(runif(100), para = c(10.2, 10.2 * 0.52, 0.52 * 3))
    
    # y <- z - x
  • 在下用过,曾搜索不少参考文献寻求p3分布拟合的实现,也曾动手写过但是无果……GAMLSS+p3分布拟合国内熊立华教授课题组使用得比较多,但是没给出具体使用方法。最后盲猜你可能是水文人吧……

    • 我在成功安装包 GPareto 之后,无法使用 library()函数加载:

      > library(GPareto)
      载入需要的程辑包:DiceKriging
      载入需要的程辑包:emoa
      Error: package or namespace load failed for ‘GPareto’ in readRDS(nsInfoFilePath):
       unknown input format
      
      > sessionInfo()
      R version 4.1.0 (2021-05-18)
      Platform: x86_64-w64-mingw32/x64 (64-bit)
      Running under: Windows 10 x64 (build 19042)
      
      Matrix products: default
      
      locale:
      [1] LC_COLLATE=Chinese (Simplified)_China.936  LC_CTYPE=Chinese (Simplified)_China.936    LC_MONETARY=Chinese (Simplified)_China.936
      [4] LC_NUMERIC=C                               LC_TIME=Chinese (Simplified)_China.936    
      
      attached base packages:
      [1] stats     graphics  grDevices utils     datasets  methods   base     
      
      other attached packages:
      [1] emoa_0.5-0.1      DiceKriging_1.6.0 mco_1.15.6       
      
      loaded via a namespace (and not attached):
      [1] compiler_4.1.0       parallel_4.1.0       tools_4.1.0          scatterplot3d_0.3-41 fortunes_1.5-4

      有人推荐了SO上的解决方案:https://stackoverflow.com/questions/6473831/readrdsfile-in-r

      看得懂的试过了没有用、还有一些方式看不懂,比如:

      Run find /usr/local/lib/R/site-library/ /usr/lib/R/library/ /usr/lib/R/site-library/ ~/.local/lib/ -iname '*rds' -a -size 0 and then delete the files found.

      也重新手动删除与加载此包,均无法解决问题,请大家支支招!

    • Liechi 谢谢!正是如此!

      你的方式也是分步解决的思路,我之前尝试用 data.table 解决,只会拆开一列,而且新列的变量难以对应:

      library(data.table)
      melt(as.data.table(tbl), 1, measure = patterns('^ak_', '^djk_'))
      • 有一个数据框:

        tbl <- structure(list(t = c(1, 2, 3, 4, 5, 6), ak_b = c(909.021709977295, 
                                                         957.828244741177, 1037.13886373248, 1122.55029956928, 1275.07072070641, 
                                                         1415.38950815257), djk_b = c(2824.00501596013, 2975.62944634725, 
                                                                                      3222.01914572632, 3487.36189890379, 3961.18824386354, 4397.10848122651
                                                         ), ak_a = c(909.021709977295, 957.828244741177, 1037.13886373248, 
                                                                     1122.55029956928, 1275.07072070641, 1415.38950815257), djk_a = c(2824.00501596013, 
                                                                                                                                      2975.62944634725, 3222.01914572632, 3487.36189890379, 3961.18824386354, 
                                                                                                                                      4397.10848122651)), row.names = c(NA, -6L), class = "data.frame")

        想把它变为如下形式:

        t   value  class1 class2
        1  909       ak        b
        2  2824     djk       b
        3  909        ak       a
        ...

        即,将列名根据_拆分之后形成的字符作为class,并对应至相应的value上,请问大家有没有方便的解法呢?

        • dupont 最好能提供一个最小可重复示例,以便具体分析。

        • 琢磨出来答案了,不要全部都用color设置,用不同的 aesthetics 即可:

          ggplot() +
            geom_point(aes(x = year, y = pop, shape = factor('observations')), data = pop, inherit.aes = TRUE) +
            geom_line(aes(x = year, y = fit, color = factor('predictions')), data = pop_pred, inherit.aes = TRUE) +
            geom_ribbon(aes(x = year, ymin = lwr, ymax = upr, fill = "0.95 confidence interval"), alpha = 0.3,
                        data = pop_pred, inherit.aes = TRUE) +
            theme_bw() +
            scale_fill_manual(values = c("0.95 confidence interval" = "grey12")) +
            scale_colour_manual(values = c("predictions" = "black")) + 
            scale_shape_manual(values = c('observations' = 16)) +
            theme(legend.title = element_blank(), text = element_text(size = 16, family = "Times New Roman"), 
                  legend.position = c(0.95, 0.05), legend.justification = c(1, 0), legend.key.height = unit(0, "pt"), 
                  legend.background = element_rect(fill = "transparent", colour = "transparent"),
                  legend.margin = margin(t = -10, r = 0, b = 0, l = 0, unit = "pt"))
        • 问题描述

          有两个数据框,pop表示实测人口,pop_pred表示模型拟合及预测的数据:其中fit为预测,lwr为 0.95 置信区间下界,upr为 0.95 置信区间上界,我想绘制一幅图,使

          • 点,表实测;
          • 线,表预测;
          • 带,表预测的置信区间

          的图例分开,应该怎么做呢?

          我的代码、数据和运行结果

          pop = structure(list(year = 1951:2020, pop = c(256.5967, 271.41, 257.12, 
          306.48, 341.55, 376.62, 411.69, 455.04, 515.03, 540.36, 486.19, 
          459.66, 470.33, 475.85, 485.06, 501.49, 505.67, 509.67, 505.49, 
          515.2, 571.55, 578.27, 579.39, 600.87, 626.22, 661.25, 679.67, 
          690.23, 747.92, 786.49, 815.6, 848.76, 873.91, 877.41, 1120.5, 
          1187.23, 1288.61, 1389.99, 1491.37, 1551.51, 1433.06, 1637.68, 
          1731.66, 1604.13, 1800.89, 1965.4, 1834.6, 1884.41, 1990.17, 
          2412.2, 2308.5, 2348.2, 2387.7, 2427.3, 2466.7, 2493.5, 2524.7, 
          2581.4, 2631, 2847, 2984, 3092, 3161, 3238, 3327, 3419, 3500, 
          3568, 3615, 3646)), row.names = c(NA, -70L), class = c("data.table", 
          "data.frame"), .internal.selfref = <pointer: 0x00000141ec3f1ef0>)
          
          pop_pred = structure(list(fit = c(215.548275690651, 227.11362269573, 239.274332172683, 
          252.058292180087, 265.494395470503, 279.612541752631, 294.443635378561, 
          310.019577881263, 326.373254749198, 343.53851578764, 361.550148380501, 
          380.44384293323, 400.256149747454, 421.024426552788, 442.786775901731, 
          465.581971621339, 489.449373511854, 514.428829489456, 540.560564389497, 
          567.885054679996, 596.442888384688, 626.274609582684, 657.420546939776, 
          689.920625836576, 723.814163792936, 759.139649047993, 795.934502342192, 
          834.234822162562, 874.075113955953, 915.488004086527, 958.503939612703, 
          1003.15087528311, 1049.45394949801, 1097.4351513484, 1147.11298122393, 
          1198.5021078675, 1251.61302513995, 1306.45171213485, 1363.01930064054, 
          1421.31175427292, 1481.31956388647, 1543.02746410021, 1606.41417593649, 
          1671.45218065217, 1738.1075298315, 1806.33969669761, 1876.10147337665, 
          1947.33891850876, 2019.99135914005, 2093.99145024978, 2169.26529457087, 
          2245.7326245579, 2323.30704745683, 2401.89635345163, 2481.40288582456, 
          2561.72397099323, 2642.75240520497, 2724.37699360587, 2806.48313638758, 
          2888.9534557789, 2971.66845681874, 3054.50721414836, 3137.34807651519, 
          3220.06938030545, 3302.55016323086, 3384.67086929108, 3466.31403631957, 
          3547.3649577893, 3627.71231109441, 3707.24874521834, 3785.87142152441, 
          3863.48250233841, 3939.98958300434, 4015.30606415666, 4089.35146203469, 
          4162.05165573932, 4233.33907137203, 4303.15280397797, 4371.43867911761, 
          4438.14925669892, 4503.24378040238, 4566.68807661461, 4628.45440725121, 
          4688.52128119289, 4746.87322928597, 4803.50054797469, 4858.3990166457, 
          4911.56959368733, 4963.01809610669, 5012.75486732131, 5060.79443745926, 
          5107.15518017678, 5151.85896964615, 5194.93084099058, 5236.39865705714, 
          5276.29278403234, 5314.64577802535, 5351.49208437785, 5386.86775111151, 
          5420.81015759915, 5453.35775924593, 5484.54984869394, 5514.42633381957, 
          5543.02753257693, 5570.39398455272, 5596.56627893693, 5621.58489847903, 
          5645.49007888787, 5668.32168304522, 5690.11908933416, 5710.92109333363, 
          5730.76582209622, 5749.6906602071, 5767.73218681403, 5784.92612282212, 
          5801.30728745845, 5816.9095634315, 5831.76586993477, 5845.90814277435, 
          5859.36732093287, 5872.1733389182, 5884.35512428248, 5895.94059973528, 
          5906.95668931319, 5917.42932810628, 5927.38347507968, 5936.84312856478, 
          5945.83134403021, 5954.37025377625, 5962.48108822884, 5970.18419853955, 
          5977.49908022674, 5984.44439762004, 5991.03800889517, 5997.29699150948, 
          6003.23766787029, 6008.87563108762, 6014.22577068133, 6019.30229812942, 
          6024.1187721591, 6028.68812369655, 6033.0226804033, 6037.13419073881, 
          6041.03384749888, 6044.73231078862, 6048.23973039683, 6051.56576754613, 
          6054.71961599927, 6057.71002250812, 6060.54530659653), lwr = c(6.23421115200421, 
          17.7093056318627, 29.7800141509885, 42.4747292642425, 55.822888135285, 
          69.8549720046533, 84.6024998994795, 100.098015867512, 116.375068983906, 
          133.468185351318, 151.412831294092, 170.245366937859, 190.002989369035, 
          210.723664586984, 232.446047497589, 255.209389253343, 279.053431324291, 
          304.018285788725, 330.144301464208, 357.471915659783, 386.041491519311, 
          415.893141143159, 447.066534918311, 479.600697751243, 513.533793176787, 
          548.902896600512, 585.743759209678, 624.090564343679, 663.975678331253, 
          705.42939795898, 748.479696812512, 793.151972708087, 839.468798288389, 
          887.449676579822, 937.110802890831, 988.464833875926, 1041.52066391261, 
          1096.28320816718, 1152.75319090306, 1210.92693676701, 1270.79616203907, 
          1332.34776221966, 1395.56359191875, 1460.42023286453, 1526.88874600417, 
          1594.93440414818, 1664.51640241403, 1735.58754483458, 1808.09390688427, 
          1881.97447531138, 1957.16076853284, 2033.57644296444, 2111.13689308126, 
          2189.74885584362, 2269.31003354868, 2349.70875338779, 2430.82368722223, 
          2512.52366150123, 2594.66759483956, 2677.10460924145, 2759.67436948568, 
          2842.20771223576, 2924.52762959247, 3006.45066780131, 3087.78878694391, 
          3168.35169842232, 3247.9496525739, 3326.39659121689, 3403.51351674831, 
          3479.13187294951, 3553.09669779354, 3625.26930814137, 3695.52931578637, 
          3763.77584864817, 3829.92794468301, 3893.92417866838, 3955.72165457131, 
          4015.29453738994, 4072.63230655782, 4127.73789464671, 4180.62584073084, 
          4231.3205479949, 4279.85469802343, 4326.26784426674, 4370.60518592644, 
          4412.91651034642, 4453.25528529316, 4491.67788035853, 4528.2428974065, 
          4563.01059218643, 4596.04237203441, 4627.40035741052, 4657.1469975683, 
          4685.34473278917, 4712.05569731511, 4737.3414584112, 4761.26278794713, 
          4783.87946356964, 4805.25009701292, 4825.4319874139, 4844.48099771288, 
          4862.45145235947, 4879.39605463891, 4895.36582200133, 4910.41003783165, 
          4924.57621814852, 4937.91009177221, 4950.4555925579, 4962.25486235202, 
          4973.34826339639, 4983.77439897723, 4993.57014119165, 5002.77066478313, 
          5011.40948607732, 5019.51850612957, 5027.12805727465, 5034.26695234644, 
          5040.96253590958, 5047.24073691616, 5053.12612226776, 5058.64195082596, 
          5063.81022747329, 5068.65175688043, 5073.18619668537, 5077.43210983518, 
          5081.40701588232, 5085.127441064, 5088.6089670262, 5091.86627808333, 
          5094.91320693027, 5097.76277874648, 5100.42725365158, 5102.91816748918, 
          5105.24637093039, 5107.42206690137, 5109.45484634954, 5111.35372237264, 
          5113.12716274165, 5114.78312085502, 5116.32906516636, 5117.77200713145, 
          5119.11852772346, 5120.37480256716, 5121.54662574431, 5122.63943232323, 
          5123.65831966586, 5124.60806756527, 5125.49315726636, 5126.31778942141, 
          5127.08590103132), upr = c(424.862340229297, 436.517939759598, 
          448.768650194378, 461.641855095932, 475.165902805721, 489.370111500608, 
          504.284770857643, 519.941139895013, 536.37144051449, 553.608846223962, 
          571.68746546691, 590.642318928601, 610.509310125873, 631.325188518592, 
          653.127504305873, 675.954553989336, 699.845315699418, 724.839373190187, 
          750.976827314785, 778.29819370021, 806.844285250066, 836.65607802221, 
          867.774558961241, 900.24055392191, 934.094534409084, 969.376401495474, 
          1006.1252454747, 1044.37907998144, 1084.17454958065, 1125.54661021407, 
          1168.52818241289, 1213.14977785813, 1259.43910070764, 1307.42062611698, 
          1357.11515955702, 1408.53938185907, 1461.70538636729, 1516.62021610251, 
          1573.28541037802, 1631.69657177883, 1691.84296573387, 1753.70716598076, 
          1817.26475995423, 1882.4841284398, 1949.32631365883, 2017.74498924704, 
          2087.68654433927, 2159.09029218294, 2231.88881139583, 2306.00842518818, 
          2381.3698206089, 2457.88880615136, 2535.47720183239, 2614.04385105964, 
          2693.49573810043, 2773.73918859866, 2854.68112318771, 2936.23032571051, 
          3018.2986779356, 3100.80230231635, 3183.6625441518, 3266.80671606097, 
          3350.16852343791, 3433.68809280959, 3517.31153951781, 3600.99004015983, 
          3684.67842006524, 3768.3333243617, 3851.91110544051, 3935.36561748717, 
          4018.64614525528, 4101.69569653545, 4184.44985022232, 4266.83627966515, 
          4348.77497938637, 4430.17913281026, 4510.95648817275, 4591.01107056601, 
          4670.2450516774, 4748.56061875113, 4825.86172007392, 4902.05560523431, 
          4977.054116479, 5050.77471811904, 5123.14127264551, 5194.08458560296, 
          5263.54274799825, 5331.46130701612, 5397.79329480688, 5462.4991424562, 
          5525.54650288411, 5586.91000294303, 5646.57094172399, 5704.51694919199, 
          5760.74161679918, 5815.24410965347, 5868.02876810357, 5919.10470518607, 
          5968.4854052101, 6016.1883277844, 6062.23452077899, 6106.64824502841, 
          6149.45661300022, 6190.68924315254, 6230.37793127379, 6268.55633972535, 
          6305.25970518585, 6340.52456521783, 6374.38850373841, 6406.88991527194, 
          6438.06778769003, 6467.9615030008, 6496.61065563107, 6524.05488755074, 
          6550.33373951466, 6575.48651764225, 6599.55217451656, 6622.56920395996, 
          6644.57554863253, 6665.60851959797, 6685.70472701043, 6704.90002109167, 
          6723.22944259013, 6740.72718194101, 6757.42654637738, 6773.35993427703, 
          6788.55881606556, 6803.05372103421, 6816.87422946917, 6830.04896952742, 
          6842.60561833262, 6854.5709068019, 6865.97062775091, 6876.82964685994, 
          6887.17191611759, 6897.02048939104, 6906.39753980259, 6915.32437862102, 
          6923.82147540381, 6931.90847915184, 6939.60424026166, 6946.92683308314, 
          6953.89357891045, 6960.52106925346, 6966.82518925401, 6972.8211411278, 
          6978.52346752698, 6983.94607473218, 6989.10225559484, 6994.00471216175
          ), year = 1951:2100), row.names = c(NA, -150L), class = c("data.table", 
          "data.frame"), .internal.selfref = <pointer: 0x00000141ec3f1ef0>)
          
          ggplot() +
            geom_point(aes(x = year, y = pop, color = 'observations'), data = pop) +
            geom_line(aes(x = year, y = fit, color = 'fit'), data = pop_pred) +
            geom_ribbon(aes(x = year, ymin = lwr, ymax = upr, fill = "0.95 confidence interval"), alpha = 0.3, data = pop_pred) +
            theme_bw() +
            scale_fill_manual(name = "", values = c("0.95 confidence interval" = "grey12")) +
            scale_colour_manual(name = "", values = c("fit" = "black", "observations" = "red"))

          预期

          我想要图例达到类似这种效果:

          • 点、线和带分开设置,图例共三行;
          • 因为带的颜色是用fill控制的,因此图例中colorfill会分开,中间隔一个空行,我想去掉这个空行,让三行如下图般紧凑,

          请问可以做到吗?

          图例.jpg

        • yihui tctcab 谢谢你们,我自己的想法有些不太规范,对 R 了解得也不太深入,所以会有这些奇奇怪怪的问题,以后注意了!

        • yihui 嗯嗯……自己的想法是数据格式在 R 中都已经处理好了,如果想要在 Rcpp 中使用能不能直接一点。可能想法很简单粗暴但是实际对于一个语言来说如果这样设计是不合理的,tctcab 既然这样就填上土吧!

          • 举个例子:

            #include <Rcpp.h>
            using namespace Rcpp;
            
            // [[Rcpp::export]]
            NumericVector f() {
              Rcpp::Environment global = Rcpp::Environment::global_env();//获取全局环境并赋值给Environment型变量global
              Rcpp::NumericVector x1 = global["x1"];
              Rcpp::NumericVector x2 = global["x2"];
              NumericVector y = x1 + x2;
              return y;
            }
            
            /*** R
            x1 <- 1:6
            x2 <- 2:7
            f()
            */

            如果 R 的环境中有很多全局变量 X1, ..., Xn,想把它们依次相加(这里只是举个例子,当然如果想做到这种相加还有别的更好的方法),而且他们都不作为输入,想只是作为全局变量,是否有比下述更简便的方法?

              Rcpp::NumericVector x1 = global["x1"];
              Rcpp::NumericVector x2 = global["x2"];
              ...
              Rcpp::NumericVector x1 = global["xn"];
            • 如题,在 R 中如果我有一个已知序列 A,想要根据不同的 R2R^2 ,以 A 为基础生成一个与 A 的相对误差服从正态分布的模拟序列 B,有什么方法吗?

              若 A 为 yty_t,B 为 yty_t',满足:
              yt=yt+ϵtyt y_t' = y_t + \epsilon_t y_t
              ϵtN(0,σ2) \epsilon_t {\sim} N(0, \sigma^2)

              • 问题

                有一组数据:

                data <- structure(list(longten = c(37.45, 42, 1300.6, 288.8, 157.05, 
                                                   487.4), tollten = c(0, 211.45, 1247.2, 0, 0, 798.4), equipten = c(0, 
                                                                                                                     0, 0, 0, 0, 0), cardten = c(110, 125, 2150, 0, 0, 570), wireten = c(0, 
                                                                                                                                                                                         380.35, 0, 0, 0, 0), longmon = c(3.7, 4.4, 18.15, 9.45, 6.3, 
                                                                                                                                                                                                                          11.8), tollmon = c(0, 20.75, 18, 0, 0, 19.25), equipmon = c(0, 
                                                                                                                                                                                                                                                                                      0, 0, 0, 0, 0), cardmon = c(7.5, 15.25, 30.25, 0, 0, 13.5), wiremon = c(0, 
                                                                                                                                                                                                                                                                                                                                                              35.7, 0, 0, 0, 0)), row.names = c(NA, -6L), class = c("tbl_df", 
                                                                                                                                                                                                                                                                                                                                                                                                                    "tbl", "data.frame"))

                它们是由*mon*ten两两对应的,想把这个数据框变形为:

                | ten | mon | type |
                | :-----: | :----: | :----: |
                | 单元格 | 单元格 | long |
                | 单元格 | 单元格 | toll |
                ...

                笨办法

                这是我用的笨办法,分别依赖于操作方式的一致才使得最后的cbind()可以对应,但是我感觉这个方法不太好。

                data_mon <- select(data, contains('mon')) %>% 
                  pivot_longer(cols = contains('mon'),
                               names_to = 'type', names_pattern = '(.*)mon',
                               values_to = 'mon')
                
                data_ten <- select(data, contains('ten')) %>% 
                  pivot_longer(cols = contains('ten'),
                               names_to = 'type', names_pattern = '(.*)ten',
                               values_to = 'ten')
                
                mon_ten <- cbind(data_mon, ten = data_ten$ten)
                
                > mon_ten
                    type   mon     ten
                1   long  3.70   37.45
                2   toll  0.00    0.00
                3  equip  0.00    0.00
                4   card  7.50  110.00
                5   wire  0.00    0.00
                6   long  4.40   42.00
                7   toll 20.75  211.45
                8  equip  0.00    0.00
                9   card 15.25  125.00
                10  wire 35.70  380.35
                11  long 18.15 1300.60
                12  toll 18.00 1247.20
                13 equip  0.00    0.00
                14  card 30.25 2150.00
                15  wire  0.00    0.00
                16  long  9.45  288.80
                17  toll  0.00    0.00
                18 equip  0.00    0.00
                19  card  0.00    0.00
                20  wire  0.00    0.00
                21  long  6.30  157.05
                22  toll  0.00    0.00
                23 equip  0.00    0.00
                24  card  0.00    0.00
                25  wire  0.00    0.00
                26  long 11.80  487.40
                27  toll 19.25  798.40
                28 equip  0.00    0.00
                29  card 13.50  570.00
                30  wire  0.00    0.00

                请问大家有无更加贴切的方法吗?

              • shang00122

                df <- data.frame(a1 = c(2.1, NA, 3.1, NA, NA, 5.1), a2 = c(1, 3, 5, 6, 7, 9))
                df$a1[is.na(df$a1)] <- df$a2[is.na(df$a1)]
                df
                #>    a1 a2
                #> 1 2.1  1
                #> 2 3.0  3
                #> 3 3.1  5
                #> 4 6.0  6
                #> 5 7.0  7
                #> 6 5.1  9

                <sup>Created on 2020-11-28 by the reprex package (v0.3.0)</sup>

              • dapengde 说出了我刚来时的心路历程,不过当时我刚好初学 Markdown,想练练手,所以并不排斥。