/DiscreteDeletor.R

http://github.com/sbotond/phylosim · R · 642 lines · 235 code · 33 blank · 374 comment · 35 complexity · d6994f94aea6d6b6ab04641e545b1e48 MD5 · raw file

  1. ##
  2. ## Copyright 2009 Botond Sipos
  3. ## See the package description for licensing information.
  4. ##
  5. ##########################################################################/**
  6. #
  7. # @RdocClass DiscreteDeletor
  8. #
  9. # @title "The DiscreteDeletor class"
  10. #
  11. # \description{
  12. # This class implements a process which performs deletions with
  13. # lengths sampled from a user-specified discrete distribution.
  14. # See \code{GeneralDeletor} for how the deletion process works.
  15. #
  16. # @classhierarchy
  17. # }
  18. #
  19. # @synopsis
  20. #
  21. # \arguments{
  22. # \item{name}{The name of the object.}
  23. # \item{rate}{The general rate.}
  24. # \item{sizes}{The deletion sizes to propose.}
  25. # \item{probs}{A vector with the probabilites of the deletion sizes.}
  26. # \item{...}{Additional arguments.}
  27. # }
  28. #
  29. # \section{Fields and Methods}{
  30. # @allmethods
  31. # }
  32. #
  33. # \examples{
  34. # # create a DiscreteDeletor process
  35. # d<-DiscreteDeletor(
  36. # name="M.D.",
  37. # rate=0.25,
  38. # sizes=c(1,2),
  39. # probs=c(1/2,1/2)
  40. # )
  41. # # get object summary
  42. # summary(d)
  43. # # set/get deletions sizes
  44. # d$sizes<-1:3
  45. # d$sizes
  46. # # set/get length probabilities
  47. # d$probs<-c(3,2,1)/6
  48. # d$probs
  49. # # plot length distribution
  50. # plot(d)
  51. #
  52. # # The following code illustrates how to use
  53. # # a DiscreteDeletor process in a simulation
  54. #
  55. # # create a sequence object, attach process d
  56. # s<-NucleotideSequence(string="AAAAAAAAAAGGGGAAAAAAAAAA",processes=list(list(d)))
  57. # # set the deletion tolerance to zero in the range 11:15
  58. # # creating a region rejecting all deletions
  59. # setDeletionTolerance(s,d,0,11:15)
  60. # # get deletion tolerances
  61. # getDeletionTolerance(s,d)
  62. # # create a simulation object
  63. # sim<-PhyloSim(root.seq=s,phylo=rcoal(2))
  64. # # simulate
  65. # Simulate(sim)
  66. # # print resulting alignment
  67. # sim$alignment
  68. # }
  69. #
  70. # @author
  71. #
  72. # \seealso{
  73. # GeneralDeletor ContinuousDeletor GeneralInDel
  74. # }
  75. #
  76. #*/###########################################################################
  77. setConstructorS3(
  78. "DiscreteDeletor",
  79. function(
  80. name="Anonymous",
  81. rate=NA,
  82. sizes=NA,
  83. probs=NA,
  84. ...
  85. ) {
  86. this<-GeneralDeletor(
  87. name=NA,
  88. rate=rate,
  89. propose.by=NA,
  90. accept.by=NA,
  91. ...
  92. );
  93. this<-extend(
  94. this,
  95. "DiscreteDeletor",
  96. .sizes=NA,
  97. .probs=NA
  98. );
  99. # Using virtual field to clear Id cache:
  100. this$name<-name;
  101. STATIC<-TRUE;
  102. if(!missing(sizes)) {
  103. this$sizes<-sizes;
  104. STATIC<-FALSE;
  105. }
  106. if(!missing(probs)) {
  107. this$probs<-probs;
  108. STATIC<-FALSE;
  109. }
  110. this$proposeBy<-function(process=NA,...){
  111. if(!exists(x="PSIM_FAST")){
  112. if( !is.numeric(process$.sizes) | !is.numeric(process$.probs) ){
  113. throw("Cannot propose deletion length because the length distribution is not defined properly!\n");
  114. }
  115. }
  116. if(length(process$.sizes) == 1){
  117. return(process$.sizes[[1]]);
  118. } else {
  119. return(sample(x=process$.sizes,size=1,replace=FALSE,prob=process$.probs));
  120. }
  121. }
  122. return(this);
  123. },
  124. enforceRCC=TRUE
  125. );
  126. ##
  127. ## Method: checkConsistency
  128. ##
  129. ###########################################################################/**
  130. #
  131. # @RdocMethod checkConsistency
  132. #
  133. # @title "Check object consistency"
  134. #
  135. # \description{
  136. # @get "title".
  137. # }
  138. #
  139. # @synopsis
  140. #
  141. # \arguments{
  142. # \item{this}{An object.}
  143. # \item{...}{Not used.}
  144. # }
  145. #
  146. #
  147. # \value{
  148. # Returns an invisible TRUE if no inconsistencies found in the object, throws
  149. # an error otherwise.
  150. # }
  151. #
  152. # @author
  153. #
  154. # \seealso{
  155. # @seeclass
  156. # }
  157. #
  158. #*/###########################################################################
  159. setMethodS3(
  160. "checkConsistency",
  161. class="DiscreteDeletor",
  162. function(
  163. this,
  164. ...
  165. ){
  166. wp<-this$writeProtected;
  167. if (wp) {
  168. this$writeProtected<-FALSE;
  169. }
  170. may.fail<-function(this) {
  171. if (is.numeric(this$sizes)) {
  172. this$sizes<-this$sizes;
  173. }
  174. else if (!is.na(this$sizes)){
  175. throw("Deletion size vector is invalid!\n");
  176. }
  177. if (is.numeric(this$probs)) {
  178. this$probs<-this$probs;
  179. }
  180. else if (!is.na(this$probs)){
  181. throw("Deletion size probability vector is invalid!\n");
  182. }
  183. }
  184. tryCatch(may.fail(this),finally=this$writeProtected<-wp);
  185. NextMethod();
  186. },
  187. private=FALSE,
  188. protected=FALSE,
  189. overwrite=FALSE,
  190. conflict="warning",
  191. validators=getOption("R.methodsS3:validators:setMethodS3")
  192. );
  193. ##
  194. ## Method: getSizes
  195. ##
  196. ###########################################################################/**
  197. #
  198. # @RdocMethod getSizes
  199. #
  200. # @title "Get the sizes of the proposed deletions"
  201. #
  202. # \description{
  203. # @get "title".
  204. # }
  205. #
  206. # @synopsis
  207. #
  208. # \arguments{
  209. # \item{this}{A DiscreteDeletor object.}
  210. # \item{...}{Not used.}
  211. # }
  212. #
  213. # \value{
  214. # A vector of integers.
  215. # }
  216. #
  217. # \examples{
  218. # # create a DiscreteDeletor object
  219. # d<-DiscreteDeletor(rate=1)
  220. # # set deletion sizes
  221. # setSizes(d,c(1,2,3))
  222. # # get deletion sizes
  223. # getSizes(d)
  224. # # set/get sizes via virtual field
  225. # d$sizes<-1:10
  226. # d$sizes
  227. # }
  228. #
  229. # @author
  230. #
  231. # \seealso{
  232. # @seeclass
  233. # }
  234. #
  235. #*/###########################################################################
  236. setMethodS3(
  237. "getSizes",
  238. class="DiscreteDeletor",
  239. function(
  240. this,
  241. ...
  242. ){
  243. this$.sizes;
  244. },
  245. private=FALSE,
  246. protected=FALSE,
  247. overwrite=FALSE,
  248. conflict="warning",
  249. validators=getOption("R.methodsS3:validators:setMethodS3")
  250. );
  251. ##
  252. ## Method: setSizes
  253. ##
  254. ###########################################################################/**
  255. #
  256. # @RdocMethod setSizes
  257. #
  258. # @title "Set the sizes of the proposed deletions"
  259. #
  260. # \description{
  261. # @get "title".
  262. #
  263. # The provided numeric vector is rounded.
  264. # }
  265. #
  266. # @synopsis
  267. #
  268. # \arguments{
  269. # \item{this}{A DiscreteDeletor object.}
  270. # \item{value}{A numeric vector.}
  271. # \item{...}{Not used.}
  272. # }
  273. #
  274. # \value{
  275. # A vector of integers (invisible).
  276. # }
  277. #
  278. # \examples{
  279. # # create a DiscreteDeletor object
  280. # d<-DiscreteDeletor(rate=1)
  281. # # set deletion sizes
  282. # setSizes(d,c(1,2,3))
  283. # # get deletion sizes
  284. # getSizes(d)
  285. # # set/get sizes via virtual field
  286. # d$sizes<-1:10
  287. # d$sizes
  288. # }
  289. #
  290. # @author
  291. #
  292. # \seealso{
  293. # @seeclass
  294. # }
  295. #
  296. #*/###########################################################################
  297. setMethodS3(
  298. "setSizes",
  299. class="DiscreteDeletor",
  300. function(
  301. this,
  302. value,
  303. ...
  304. ){
  305. .checkWriteProtection(this);
  306. if (missing(value)) {
  307. throw("No new value provided!\n");
  308. } else if (!is.numeric(value)) {
  309. throw("The new value must be numeric vector!\n");
  310. } else {
  311. if(length(value) == 0 ) {
  312. warning("Deletion size vector has zero length!\n");
  313. }
  314. this$.sizes<-round(value);
  315. }
  316. return(invisible(this$.sizes));
  317. },
  318. private=FALSE,
  319. protected=FALSE,
  320. overwrite=FALSE,
  321. conflict="warning",
  322. validators=getOption("R.methodsS3:validators:setMethodS3")
  323. );
  324. ##
  325. ## Method: getProbs
  326. ##
  327. ###########################################################################/**
  328. #
  329. # @RdocMethod getProbs
  330. #
  331. # @title "Get the deletion length probabilities"
  332. #
  333. # \description{
  334. # @get "title".
  335. # }
  336. #
  337. # @synopsis
  338. #
  339. # \arguments{
  340. # \item{this}{A DiscreteDeletor object.}
  341. # \item{...}{Not used.}
  342. # }
  343. #
  344. # \value{
  345. # A numeric vector with the deletion length probabilities.
  346. # }
  347. #
  348. # \examples{
  349. # # create a DiscreteDeletor object
  350. # d<-DiscreteDeletor(rate=1, sizes=1:3)
  351. # # set/get length probabilities
  352. # setProbs(d,c(1/3,1/3,1/3)) # equal probabilites
  353. # getProbs(d)
  354. # # set/get length probabilities via virtual field
  355. # x<-c(2,2,1)
  356. # # normalize x
  357. # x<-x/sum(x)
  358. # d$probs<-x
  359. # d$probs
  360. # }
  361. #
  362. # @author
  363. #
  364. # \seealso{
  365. # @seeclass
  366. # }
  367. #
  368. #*/###########################################################################
  369. setMethodS3(
  370. "getProbs",
  371. class="DiscreteDeletor",
  372. function(
  373. this,
  374. ...
  375. ){
  376. this$.probs;
  377. },
  378. private=FALSE,
  379. protected=FALSE,
  380. overwrite=FALSE,
  381. conflict="warning",
  382. validators=getOption("R.methodsS3:validators:setMethodS3")
  383. );
  384. ##
  385. ## Method: setProbs
  386. ##
  387. ###########################################################################/**
  388. #
  389. # @RdocMethod setProbs
  390. #
  391. # @title "Set the deletion length probabilities"
  392. #
  393. # \description{
  394. # @get "title".
  395. #
  396. # The \code{sizes} virtual field must be set before setting the length probabilities.
  397. # The length of the provided numeric vector must match with the length of the vector
  398. # stored in the \code{sizes} virtual field. The vector is rescaled if the values do not
  399. # sum to one and a warning is issued.
  400. # }
  401. #
  402. # @synopsis
  403. #
  404. # \arguments{
  405. # \item{this}{A DiscreteDeletor object.}
  406. # \item{value}{A numeric vector containg the length probabilities.}
  407. # \item{...}{Not used.}
  408. # }
  409. #
  410. # \value{
  411. # The vector of probabilities.
  412. # }
  413. #
  414. # \examples{
  415. # # create a DiscreteDeletor object
  416. # d<-DiscreteDeletor(rate=1, sizes=1:3)
  417. # # set/get length probabilities
  418. # setProbs(d,c(1/3,1/3,1/3)) # equal probabilites
  419. # getProbs(d)
  420. # # set/get length probabilities via virtual field
  421. # x<-c(2,2,1)
  422. # # normalize x
  423. # x<-x/sum(x)
  424. # d$probs<-x
  425. # d$probs
  426. # }
  427. #
  428. # @author
  429. #
  430. # \seealso{
  431. # @seeclass
  432. # }
  433. #
  434. #*/###########################################################################
  435. setMethodS3(
  436. "setProbs",
  437. class="DiscreteDeletor",
  438. function(
  439. this,
  440. value,
  441. ...
  442. ){
  443. .checkWriteProtection(this);
  444. if (missing(value)) {
  445. throw("No new value provided!\n");
  446. }
  447. else if(!is.numeric(this$.sizes)) {
  448. throw("Cannot set probabilities because indel size vector is not defined!\n");
  449. }
  450. else if (!is.numeric(value)) {
  451. throw("The new value must be a numeric vector!\n");
  452. }
  453. else if(length(value) != length(this$.sizes)) {
  454. throw("The length of the probabilities vector must be the same as the length of the deletion sizes vector");
  455. }
  456. else if( length(value[value < 0 ]) != 0 ) {
  457. throw("The elements of the probability vector must not be negative!\n");
  458. }
  459. else {
  460. if(!isTRUE(all.equal(sum(value),1.0))){
  461. value<-value/sum(value);
  462. warning("The provided values were rescaled in order to sum to one!\n");
  463. }
  464. this$.probs<-value;
  465. }
  466. },
  467. private=FALSE,
  468. protected=FALSE,
  469. overwrite=FALSE,
  470. conflict="warning",
  471. validators=getOption("R.methodsS3:validators:setMethodS3")
  472. );
  473. ##
  474. ## Method: summary
  475. ##
  476. ###########################################################################/**
  477. #
  478. # @RdocMethod summary
  479. #
  480. # @title "Summarize the properties of an object"
  481. #
  482. # \description{
  483. # @get "title".
  484. # }
  485. #
  486. # @synopsis
  487. #
  488. # \arguments{
  489. # \item{object}{An object}
  490. # \item{...}{Not used.}
  491. # }
  492. #
  493. # \value{
  494. # Returns a PSRootSummary object.
  495. # }
  496. #
  497. # \examples{
  498. # # create an object
  499. # a<-DiscreteDeletor(rate=1,sizes=c(1,2),probs=c(1/2,1/2))
  500. # # get a summary
  501. # summary(a)
  502. # }
  503. #
  504. # @author
  505. #
  506. # \seealso{
  507. # @seeclass
  508. # }
  509. #
  510. #*/###########################################################################
  511. setMethodS3(
  512. "summary",
  513. class="DiscreteDeletor",
  514. function(
  515. object,
  516. ...
  517. ){
  518. this<-object;
  519. .addSummaryNameId(this);
  520. expected.length<-NA;
  521. sd<-NA;
  522. if( is.numeric(this$sizes) & is.numeric(this$probs)) {
  523. expected.length<-weighted.mean(this$sizes, this$probs);
  524. sd<-sqrt(sum( (this$sizes - expected.length)^2 * this$probs ));
  525. }
  526. this$.summary$"Expected deletion length"<-expected.length;
  527. this$.summary$"Standard deviation of deletion length"<-format(sd);
  528. NextMethod();
  529. },
  530. private=FALSE,
  531. protected=FALSE,
  532. overwrite=FALSE,
  533. conflict="warning",
  534. validators=getOption("R.methodsS3:validators:setMethodS3")
  535. );
  536. ##
  537. ## Method: plot
  538. ##
  539. ###########################################################################/**
  540. #
  541. # @RdocMethod plot
  542. #
  543. # @title "Plot the deletion length distribution"
  544. #
  545. # \description{
  546. # @get "title".
  547. # }
  548. #
  549. # @synopsis
  550. #
  551. # \arguments{
  552. # \item{x}{A DiscreteDeletor object.}
  553. # \item{...}{Not used.}
  554. # }
  555. #
  556. # \value{
  557. # The DiscreteDeletor object (invisible).
  558. # }
  559. #
  560. # \examples{
  561. # d<-DiscreteDeletor(
  562. # name="MyDiscDel",
  563. # sizes=1:6,
  564. # probs=c(0.25000000, 0.16666667, 0.08333333, 0.08333333, 0.16666667, 0.25000000)
  565. # )
  566. # # plot the deletion length distribution
  567. # plot(d)
  568. # }
  569. #
  570. # @author
  571. #
  572. # \seealso{
  573. # @seeclass
  574. # }
  575. #
  576. #*/###########################################################################
  577. setMethodS3(
  578. "plot",
  579. class="DiscreteDeletor",
  580. function(
  581. x,
  582. ...
  583. ){
  584. this<-x;
  585. if( !is.numeric(this$sizes) | !is.numeric(this$probs) ){
  586. warning("Deletion length distribution is not defined properly! Nothing to plot here!\n");
  587. return();
  588. }
  589. plot.default(
  590. x=this$sizes,
  591. y=this$probs,
  592. col=c("blue"),
  593. lwd=2,
  594. type="h",
  595. main=paste("Deletion size distribution for:",this$id),
  596. xlab="Size",
  597. ylab="Probability",
  598. ylim=c(0,1),
  599. xaxt="n"
  600. );
  601. axis(side=1, at=this$sizes, labels=this$sizes);
  602. return(invisible(this));
  603. },
  604. private=FALSE,
  605. protected=FALSE,
  606. overwrite=FALSE,
  607. conflict="warning",
  608. validators=getOption("R.methodsS3:validators:setMethodS3")
  609. );