Question

Anyone finished with PS 12? I'm working on Problem 5. I think I got everything working properly (the resistant virus chart from problem 4 behaves a I would expect), but when I make the histogram for the various treatment schedules, I get 0% recovery rate for every timelapse except t = 0, which obviously has a 100% recovery rate. This sounds like it could be reasonable, but it doesn't seem that interesting, so I wanted to see if anyone else got these results. Please let me know!

Answer 1

i'm winding up somewhere between 500-600 final population no matter what timestep I introduce the drug

Answer 2

oh, and a similar issue with #7: when I introduce both drugs at the same time, the "recovery rate" is 100%; any other timelapse between the two drugs and I get similar results: approx 500-600 final population every time. Does this sound right?

Answer 3

Hi. I think I just finished ps12 too, but Im not getting what you mentioned here. first, I don't think #5 should be zero, because you run the trial for 150 steps before giving any medicine. during this time, there should probably be a few viruses that mutationally become immune. Are you ading the drug before the first 150? However, on #7, I think that for 0 timesteps, there should be a complete cure, because the iruses do not have time to become immune to the drugs and therefore cannot reproduce, therefore only leaving the possibility to die through the clearProb(). im not getting the 500-600 range but instead Im getting maxed out 1000 (Im unsure about this) If it's ok with you, I'd like to compare code (mine is atached)

Answer 4

hey, I thought problem 5 was supposed to have you introduce the drug at various timesteps, with one of the trials being t=0 (basically where you have 100 viruses and none of them are immune to guttagonol). So I'm not getting 0 for all of problem 5, just the one where I introduce the drug before running any simulation timesteps. My concern is that I'm not getting any interesting spread on the histograms; they're either all at 0 or all in the 500-600 range. I attached my code.

Answer 5

uhh, sorry i didn't comment my code much, just ask on here if there's anything I did and you don't understand why.

Answer 6

edit: i do get the same results as you do on #5 except I get the maxed out 1000 for everything other than 0 steps before applying the medicine. sorry I got this mixed up with another one :P forgive me

Answer 7

oh ok.. it makes me feel better that you're getting similar results. However, I think you may be missing an element if you're maxing out at 1000. There's the factor of population density where the rate of reproduction goes down as the population increases. This should be factored into the repro rate, and I believe it plays a significant enough role where by the time you get to about 600 viruses, the clear rate and virus repro rate hit an equilibrium where it fluctuates around 600. Are you accounting for that in your code? That was an element of Problem 1.

Answer 8

I think I did, but I'll recheck

Answer 9

I have it compared to clearProb such as the following:
reproduction = self.maxBirthProb * (1- popDensity)
if reproduction > self.clearProb:
return SimpleVirus(self.maxBirthProb, self.clearProb)

Answer 10

thats in my simpleVirus class. Did you compare it to random.random()?

Answer 11

I thought I was supposed to compare it to clearProb... that might be my bug

Answer 12

yeah, oh that might be an issue... I don't think you're supposed to. I think the clear rate and the repro rate are two separate tests that are supposed to be applied to each virus on each timestep. Make sense? i could be wrong, but that's how I read it, and it seems more intuitive than comparing the two rates to each other

Answer 13

hey i'm going to sleep now -- good luck, i'll be back in the morning

Answer 14

yeah It is getting late. I'll try to work on this a little more tonight. ttyl

Answer 15

So I was looking over My code last night and realized that even though I had a messed up by comparing reproduction probability to clearProb in simpleVirus, I had done the correct thing in ResistantVirus. Since that should override the reproduce in simpleVirus, I don't thinkthat's what's causing the discrepancy.

Answer 16

Ahh OK, I'll take a look and try to figure out why you're getting to 1000 viruses.

Answer 17

I think the issue is that you're creating a separate list of babies on each time step which you then append to the total viruses list. Because of this, the pop density never has the chance to buffer the population because you're taking a new list of babies (the size of which was affected by the pop density) and then just adding it to the viruses list, capping it at 1000.
I think the way to do this is to just use the virus.reproduce() function to directly add to the patient's viruses list. That way you won't be adding smaller lists to the existing list, and the pop density ratio can actually put pressure on your overall viruses list. You won't need the hard cap at 1000 either.

If you think about it, the equilibrium makes sense, because the ratio of pop density will cause the likelihood of reproduction to decrease until it's approximately the same as the clear rate (0.05), at which point they will balance out and the virus population will just fluctuate around that number.

Answer 18

aaaaaaaah I think Your right. Let me see if I can fix it

Answer 19

However, the instructions specifically state to return the offspring in reproduce, not append to the virus list so I'll have to find a way around that part

Answer 20

ok so Im on #2 because of my messed-up-ness (lol) for the portion that you mentioned I tried this:
for j in self.viruses:
populationDensity = float(len(self.viuses))/ self.maxPop
try: self.viruses.append(j.reproduce(populaionDensity)
except NochildException: pass

but I got a flat 500... line

Answer 21

so I tried this ( I know I should be using dpaste, but im running the code on a different computer, forgive me)
babies = []
for j in self.viruses:
populationDensity = float(len(self.viruses))+len(babies)/self.maxPop
try: babies.append(j.reproduce(populationDensity))
except NoChildException: pass

for baby in babies:
if len(self.viruses) < maxPop: self.viruses.append(baby)

Answer 22

Im getting the feeling that this isnt right....

Answer 23

wait a sec... order of operations...duh

Answer 24

I got my code working and It seems to be doing fine but I don't think Im getting the exactly same answers as you are... not sure. My histograms are hard to read, and I never did quite understand hist.any suggestions? oh and this probably isn't the most treamlined way to do it but it's what I've git right now

Answer 25

(>_<)

Answer 26

ok, looking at it now. Your chart that comes up for problem 4 concerns me.. I believe that was is supposed to happen with the resistant population is that it is supposed to increase exponentially from the beginning of the trial, without being affected by the introduction of guttagonol. Your resistant population seems to decrease dramatically after the introduction of guttagonol, which is something I have to look into.

So basically the resistant population, if you think about it, will start out as a small portion of the total population and grow very slowly, but then once guttagonol is introduced, the non-resistant ones quickly die off and the resistant ones continue to reproduce, so very shortly after the drug is administered you have a population that is close to 100% resistant. I'll look through it again and try to figure out why your resistant viruses are dying off after introducing guttagonol.

Also, I noticed your mutation rate on problem 4 was 0.45 -- were you just playing around with some different values there?

I think most of your histograms are now more or less similar to mine, i.e. pretty much all your results come in around 500-600 population. So you still have the issue that originally concerned me: the lack of variety. Like I thought maybe if you introduced the drug after t=75 for example and ran the test a bunch of times, you'd have a chance of actually killing off the population on some trials, but it seems like no matter when you introduce the drug (other than t=0, when there are no resistant viruses), you wind up with the population around 500-600. I'd recommend removing the range 0-1000 parameter on your histograms -- it will automatically set up the bins to give you a little more detailed view. Other than that your histograms look fine to me...

Answer 27

oh, I'm gonna try to figure out why your resistant population is going down a little later this afternoon... i gotta go eat smth now

Answer 28

ok thanks I'll probably go to the y and swim then take a new look at this

Answer 29

and yes I was playing with the numbers a bit

Answer 30

P.S. (no not as in playstation)
When I changed the mutation raet back to .005 it seems to do the jump you mentioned (which is logically correct). As for your concern about the virus not dying out when timesteps > 75 on #5, I was thinking that since you give the virus an opportunity to mutate and become resistant to the virus, you can't exactly kill it all off. Just m intuition.

Also, looking at your code, I noticed that your code for the problems that need to be spitting out multiple graphs/charts. The pdf doesn't seem to clarify this but when you reed the given comments under each problem in the code, I inferred that you were supposed to display all the charts at once, since it says 'Histograms' (plural). Just a side note :)

oh and if your just going to do one chart, I don't think you will need pylab.figure() since you don't need to worry about everything going onto the same chart/graph

yeah, long post script

Answer 31

yeah, OK I just checked it out too... that mutation was the only thing wrong with the resistant viruses. Thanks for the tips on the chart display. Yeah I was just displaying one chart and then changing the variables to see the different ones. OK! Now off for some playstation...

Answer 32

Thanks for all your help, and I hope I wasn't a nuisance. Now Im taking harvard's opencourseware cs50 but Im going to have to cram it into like 2 weeks because college is going to start next month :) Don't worry though, I will do the final exam.

Answer 33

sounds good! good luck with your studies. Haven't see the harvard OCW... i"ll check that out myself. I moved on to 6.001 here, but it's a lecture from 1985, so I'm trying to convince myself it's still relevant. I think it focuses on technology-independent computing concepts, so it looks pretty good.

Answer 34

I downloaded all the 6.001 stuff, but erased it because of the oldschoolness and I don't think I exactly need scheme (plus I couldn't get scheme to run correctly in the first place haha). Some sources ( http://www.cplusplus.com/forum/lounge/13001/) say the course has been replaced at MIT, so Im not sure if it's worth the time investment, but that's personal bias I guess. Oh and In parting, good luck to you also!